Diffuse-Charge Dynamics in Electrochemical Systems 
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The response of a model micro-electrochemical system to a time-dependent applied voltage is 
analyzed. The article begins with a fresh historical review including electrochemistry, colloidal 
science, and microfluidics. The model problem consists of a symmetric binary electrolyte between 
parallel-plate, blocking electrodes which suddenly apply a voltage. Compact Stern layers on the 
electrodes are also taken into account. The Nernst-Planck-Poisson equations are first linearized 
and solved by Laplace transforms for small voltages, and numerical solutions are obtained for large 
voltages. The "weakly nonlinear" limit of thin double layers is then analyzed by matched asymptotic 
expansions in the small parameter e = \d/L, where Ad is the screening length and L the electrode 
separation. At leading order, the system initially behaves like an RC circuit with a response time of 
\dL/D (not \\,/D), where D is the ionic diffusivity, but nonlinearity violates this common picture 
and introduce multiple time scales. The charging process slows down, and neutral-salt adsorption by 
the diffuse part of the double layer couples to bulk diffusion at the time scale, l? jD. In the "strongly 
nonlinear" regime (controlled by a dimensionless parameter resembling the Dukhin number), this 
effect produces bulk concentration gradients, and, at very large voltages, transient space charge. 
The article concludes with an overview of more general situations involving surface conduction, 
multi-component electrolytes, and Faradaic processes. 
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I. INTRODUCTION 

There is rapidly growing interest in micro- 
electrochemical or biological systems subject to time- 
dependent applied voltages or currents. For example, AC 
voltages applied at microelectrodes can be used to pump 
liquid electrolytes [llllilllllllllili[l3[llf, 
to separate or self-assemble colloidal parti- 
cles [12, 13, Q m 0,6, 17 181 and to manipulate 
biological cells and vesicles 0, 1^- Conversely, 
oscillating pressure-driven flows can be used to produce 
frequency-dependent streaming potentials to probe the 
structure of porous media [l^, IH, • 

A common feature of these diverse phenomena is the 
dynamics of diffuse charge in microscopic systems. Al- 
though the macroscopic theory of neutral electrolytes 
with quasi-equilibrium double layers is very well de- 
veloped in electrochemistry [2^ and colloidal sci- 
ence [23, 13 microscopic double-layer charging at 
subdifFusive time scales is not as well understood. Al- 
though much progress has been made in various disjoint 
communities, it is not so widely appreciated, and some 
open questions remain, especially regarding nonlinear ef- 
fects. The goals of this paper are, therefore, (i) to review 
the relevant literature and (ii) to analyze a basic model 
problem in considerable depth, highlighting some new re- 
sults and directions for further research. 

To illustrate the physics of diffuse-charge dynamics, 
consider the simplest possible case sketched in Fig. ^ a 
dilute z:z electrolyte suddenly subjected to a DC voltage, 
2V , by parallel-plate blocking electrodes separated by 2L. 
Naively, one might assume a uniform bulk electric field, 
E — V/L, but the effect of the applied voltage is not 
so trivial. Ions migrate in the bulk field and eventually 
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FIG. 1: Sketch of the model problem. A voltage 2V is sud- 
denly applied to a dilute, symmetric, binary electrolyte be- 
tween parallel-plate, blocking electrodes separated by 2L. 



screen it completely (since "blocking electrodes" do not 
support a Faradaic current). 

What is the characteristic time scale of this response? 
For charge relaxation, one usually quotes the time, td — 
Xjj/D, for diffusion with a diffusivity D across one Debye 
screening length. 
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where Cb is the average solute concentration, k Boltz- 
mann's constant, T the temperature, e is the electronic 
charge, and e the permittivity of the solvent ^7, 28, 29|. 
The Debye time, td, is a material property of the elec- 
trolyte, which for aqueous solutions (Ad « 1 — 100 nm. 
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D K, IQ^ /im^/s) is rather small, in the range of ns to 
/iS. More generally, when Faradaic reactions occur (for 
a non-blocking electrode), diffuse charge may also vary 
on the much slower, geometry-dependent scale for bulk 
diffusion given by tl = l? jD^ proportional to the square 
of the electrode separation. 

These two relaxation times, ty) for the charge density 
and for the concentration, are usually presented as 
the only ones controlling the evolution of the system, 
e.g. as in the recent textbooks of Hunter |23| (Ch. 8) and 
Lyklema (Chs. 4.6c). Dimensional analysis, however, 
allows for many other time scales obtained by combining 
these two, such as the harmonic mean. 
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proportional to the electrode separation (not squared). 
Below, we will show that this is the primary time scale 
for diffuse-charge dynamics in electrochemical cells, al- 
though T£), tl, and other time scales involving surface 
properties also play important roles, especially at large 
voltages (even without Faradaic processes). The same 
applies to highly polarizable or conducting colloidal par- 
ticles, where L is the particle size. 

Although the basic chargin g ti me, Tc, is familiar in 
several scientific communities |3ll Isa . Issl Is^ , it is not 
as widely known as it should be. Recently, it has been 
rediscovered as the (inverse) frequency of "AC pumping" 
at patterned-surface micro-electrodes P, @. As in the 
past, its theoretical justification has sparked some con- 
troversy ^35, 36| related to the applicability of classical 
circuit models |33,|33 in which Tc arises as the "RC time" 
of a bulk resistor in series with a double-layer capacitor 
(see below). 

Here, we attempt to unify and modestly extend a large 
body of prior work on diffuse-charge dynamics in the con- 
text of our model problem, paying special attention to 
effects which undermine the classical circuit approxima- 
tion. Going beyond most previous mathematical stud- 
ies, we allow for compact-layer capacitance, bulk con- 
centration polarization, and large voltages outside the 
linear regime. For the nonlinear analysis, the method 
of matched asymptotic expansions [s^ ^| must be 
adapted for multiple time scales at different orders of 
the expansion, so the problem also presents an oppor- 
tunity for mathematicians to develop time-dependent 
boundary-layer theory. 

We begin in section^by reviewing some of the relevant 
literature on electrochemical relaxation. In section llTTI we 
state the mathematical problem for a suddenly applied 
DC voltage, and in section Hvl we analyze the linear re- 
sponse using Laplace transforms. In section we non- 
dimensionalize the problem and describe numerical solu- 
tions, used to test our analytical approximations. In sec- 
tion we derive uniformly valid asymptotic expansions 
in the "weakly nonlinear" limit of thin double layers and 
discuss the connection with circuit models. Apparently 
for the first time (for this problem), in section [VIII we 



analyze higher-order corrections, and in section rVIIII we 
briefly discuss the "strongly nonlinear" regime at large 
voltages, where the expansions are no longer valid. In 
section Hxl we conclude by briefly discussing extensions 
to higher dimensions, general electrolytes, and Faradaic 
processes and posing some open questions. 



II. HISTORICAL REVIEW 
A. Electrical Circuit Models 

In electrochemistry, the most common theoretical ap- 
proach is to construct an equivalent electrical circuit, 
whose parameters are fit to experimental impedance 
spectra or pulsed-voltage responses, as recently reviewed 
by Macdonald [s^ and Geddes |3l| . The basic idea of an 
equivalent circuit is apparently due to Kohlrausch 42] in 
1873, and the first mathematical theory of Kohlrausch's 
"polarization capacitance" was given by Warburg at the 
end of the nineteenth century jl^, 1131 ■ Warburg argued 
that AC electrochemical response is dominated by pure 
diffusion of the active species and can be described a bulk 
resistance in series with a frequency-dependent capaci- 
tance, which combine to form the "Warburg impedance" . 

Earlier, Helmholtz Q| had suggested that the 
solid-electrolyte interface acts like a thin capacitor, for 
which he apparently coined the term, "double layer" [23 . 
In 1903 Kriiger i8\ unified Warburg's bulk impedance 
with Helmholtz' double-layer capacitor in the first com- 
plete AC circuit model for an electrochemical cell, which 
forms the basis for the modern "Randies circuit" [i^ . 
In this context, the relaxation time for charging of the 
double layers has been known to depend on the electrode 
separation, via the bulk resistance, for at least a century. 

The study of diffuse charge in the double layer was ini- 
tiated in the same year by Gouy |0 , who suggested that 
excess ionic charge in solution near the electrode could 
be viewed as a capacitance, Co — s/Xd- He was also the 
first to derive Equation ^ for A u (obviously with differ- 
ent notation) in his original theory of the diffuse double 
layer in equilibrium [sOllSlj . With the availability of Ein- 
stein's relation for the mobility, /i — D/kT, at that 
time, the DC bulk resistance (per unit area) could have 
been calculated as 
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(for a symmetric binary electrolyte of equal mobilities), 
where J is the current density and 
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is the bulk conductivity. Therefore, basic time scale in 
Eq. ((2Jl has essentially been contained in circuit models 
since roughly 1910 as the relaxation time, 
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although Tc was not stated exphcitly as XdL/D for per- 
haps another fifty years [3ll |. 

Today, Gouy's screening length bears the name of De- 
bye, who rederived it in 1923 as part of his seminal 
work with Hiickel jsj on charge screening in bulk 
electrolytes, using an equivalent formalism. Debye and 
Hiickel solved for the spherical screening cloud around 
an ion, and, due to the low potentials involved, they lin- 
earized the transport equations, allowing them to han- 
dle general electrolytes. When Gouy [U considered the 
identical problem of screening near a flat, blocking elec- 
trode more than a decade earlier, he obtained exact so- 
lutions to full nonlinear equations for the equilibrium 
potential profile in several cases of binary electrolytes, 
z+/z- = 1, 2, and ^, where and Z- are the cation 
and anion charge numbers, res pec tively. 

A few years later. Chapman |53 independently derived 
Gouy's solution for a univalent electrolyte, z+ = z- = I, 
the special case of "Gouy-Chapman theory" for which 
they are both primarily remembered today. Chapman 
also gave a simple form for the charge-voltage relation 
of the diffuse-layer capacitor in this case, which, upon 
differentiation, yields a simple formula for the nonlinear 
differential capacitance of the diffuse layer, 

where C is the voltage across the diffuse layer in thermal 
equilibrium. (Here, we include the trivial extension to a 
general z:z electrolyte.) Combining Eqs. Q and JHJl, we 
also obtain the basic relaxation time, Tc, in Eq. ^ mul- 
tiplied by a potential-dependent factor in the usual case 
of nonzero equilibrium zeta potential (in the absence of 
an applied voltage). This factor may be neglected in the 
Debye-Hiickel limit of small potentials, C kT/ze, but 
it becomes important at large potentials and generally 
slows down the final stages of double-layer charging. 

More sophisticated models of the double layer were 
proposed by many subsequent authors .57] and in- 
corporated into AC circuit models for electrochemi- 
cal cells |2^ |5^ [S^. Naturally, the original ideas of 
Helmholtz and Gouy were eventually combined into a 
coherent whole. In 1924, Stern suggested decom- 
posing the double layer into a "compact" (Helmholtz) 
part within a molecular distance of the surface and a 
"diffuse" (Gouy) part extending into the solution at the 
scale of the screening length. Physically, the compact 
layer is intended to describe ions (at the outer Helmholtz 
plane) whose solvation molecules are in contact with the 
surface, although specifically adsorbed ions (themselves 
in contact with the surface) may also be included [6ll |. 
Regardless of the precise microscopic picture, however. 
Stern introduced the compact layer as an intrinsic sur- 
face capacitance, which cuts off the divergent capacitance 
of the diffuse layer, Eq.©, at large zeta potentials. 

Using this model of two capacitors in series and ne- 
glecting specific adsorption, Grahame '62] applied Gouy- 
Chapman theory for the diffuse part and inferred the non- 



linear differential capacitance of the compact part from 
his famous experiments on electrified liquid-mercury 
drops. Macdonald [s^l then developed a mathematical 
model for double layers at metal electrodes by viewing 
the compact layer as a parallel-plate capacitor, as we 
do below, although he also allowed its thickness and ca- 
pacitance to vary due to electrostriction and dielectric 
saturati on .641 ■ The reader is referred to various recent 
reviews |37l l38l Isgj to learn how other effects neglected 
below, such as specific adsorption and Faradaic processes 
for non-blocking electrodes, have been included empiri- 
cally in modern circuit models. 

In spite of a century of research, open questions remain 
about the applicability of circuit models [s^, and even 
the most sophisticated fits to experimental data still suf- 
fer from ambiguities .37]. One problem is the somewhat 
arbitrary distinction between the diffuse layers and the 
bulk electrolyte, which in fact comprise a single, continu- 
ous region. Even accepting this partition, it is clear that 
the non-uniform evolution of ionic concentrations in both 
regions cannot be fully captured by homogeneous circuit 
elements [65| . Another problem is the further partition- 
ing of the double layer into two (or more) poorly defined 
regions at atomic lengths scales, where macroscopic con- 
tinuum theories (e^. for dielectric response) are of ques- 
tionable validity |66| . 



B. Microscopic Transport Models 

An alternative theoretical approach, pursued below, is 
to solve the time-dependent Nernst-Planck equations [63, 
|6^ |6^ for ionic transport across the entire cell (out- 
side any molecular-scale compact layers) without distin- 
guishing between the diffuse-charge layers and the quasi- 
neutral bulk. Because this "phenomenological" ap- 
proach requires solving Poisson's equation for the mean- 
field electrostatic potential (self-consistently generated 
by the continuum charge density) down to microscopic 
(and sometimes atomic) length scales, it lacks the ther- 
modynamic justification of traditional macroscopic the- 
ories based on bulk electroncutrality and electrochem- 
ical potentials [26|. Nevertheless, it addresses time- 
dependent charge-relaxation phenomena, which do occur 
in real systems, with fewer ad hoc assumptions than cir- 
cuit models and thus may be considered closer to first 
principles. The use of the Nernst-Planck equations at 
scales smaller than the screening length (but still larger 
than atomic dimensions) is also supported by the success 
of Gouy-Chapman theory in predicting the diffuse-layer 
capacitance in a number of experimental systems (e.g. 
Refs. IH^iliH), because the theory is based on the steady- 
state Nernst-Planck equations for thermal equilibrium. 
The main difficulty in working with the Nernst-Planck 
equations, aside from mathematical complexity, is per- 
haps in formulating appropriate boundary conditions at 
the electrode surface, just outside any compact layers. 

Although the response to a suddenly applied DC volt- 



4 



age has been considered by a few authors in the hn- 
ear [t^, IS [z3 and nonhnear [s^, |23 regimes, as we 
also do below, much more analysis has been reported for 
the case of weak AC forcing, where the equations are lin- 
earized and the time dependence is assumed to be sinu- 
soidal. These simplifications are made mainly for analyt- 
ical convenience, although they have direct relevance for 
the interpretation of impedance spectra. An early anal- 
ysis of this type was due to Ferry [s^l, who considered 
the response of a semi-infinite electrolyte to an oscillat- 
ing charge density applied at a electrode surface. Ferry's 
treatment is formally equivalent to the classical theory 
of dielectric dispersion in bulk electrolytes 0, [t^ ■ Nat- 
urally, in both cases the same time scale, td = }?jj/D, 
arises, and the relaxation of the double layer has no de- 
pendence on the macroscopic geometry. 

Ferry's analysis of a single electrode is consistent with 
the common intuition that double-layer charging should 
be a purely microscopic process, but one might won- 
der how the electrode could draw charge "from infinity" 
when an infinite electrolyte has infinite resistance. In- 
deed, as emphasized by Buck ^] and Macdonald [sj 
and confirmed by detailed comparisons with experimen- 
tal impedance spectra. Ferry's analysis is fundamentally 
flawed, starting from the boundary conditions: It is not 
possible to control the microscopic charge density at an 
electrode surface and neglect its coupling to bulk trans- 
port processes; instead, one imposes a voltage relative 
to another electrode and observes the resulting current 
(or vice versa), while the surface charge density evolves 
self-consistently. 

Buck eventually corrected Ferry's analysis to ac- 
count for the missing "IR drop" across two electrodes 
which imposes the initial surface charge density (or, al- 
ternatively, voltage). Nevertheless, the physical picture 
of a double-layer responding locally to "charge injec- 
tion" , independent of bulk transport processes, persists 
to the present day. For example, recent textbooks on col- 
loidal science (Hunter p. 408; Lyklema [l^, p. 4.78) 
present a slightly different version of Ferry's analysis (at- 
tributed to O'Brien) as the canonical problem of "double- 
layer relaxation": the response of a semi- infinite elec- 
trolyte to a suddenly imposed, constant surface charge 
density. This gives some insight into high frequency di- 
electric dispersion of non-polarizable colloids (the usual 
case), but it is not relevant for polarizable particles and 
electrodes. Three decades after Buck and Macdonald, 
it is worth re-emphasizing the fundamental coupling of 
double-layer charging of bulk transport in finite, polariz- 
able systems. 

The mathematical theory of AC response for a fi- 
nite, two-electrode system began with Jaffe's analysis for 
semiconductors (I^k iZZij and was extended to liquid elec- 
trolytes by Chang and Jaffe [t^ • A number of restrictive 
assumptions in these studies, such as a uniform electric 
field, were relaxed by Macdonald [t^ for semiconduc- 
tors and electrolytes and independently by Friauf [s^ 
for ionic crystals. These authors, who gave perhaps the 



first complete mathematical solutions, also allowed for 
bulk generation/recombination reactions, which are cru- 
cial for electrons and holes in semiconductors. Subse- 
quent authors mostly neglected bulk reactions in studies 
of liquid ^ 73, 81. M, Ml and solid JiMj electrolytes, 
while focusing on other effects, such as arbitrary ionic va- 
lences and the compact layer. 

Although it is implicit in earlier work, Macdonald [sj 
first clearly identified the geometry-dependent time scale, 
Tc = ^dL/D^ (in this form) as governing the relaxation 
of an electrochemical cell. It was also derived indepen- 
dently by Kornyshev and Vorontyntsev [tiL It^ in the 
Russian literature on solid electrolytes with one mobile 
ionic species [s^- With Itskovich [s^, these authors 
also modeled the compact-layer capacitance via a mixed 
Dirichlet-Neumann condition on the Nernst-Planck equa- 
tions. This classical boundary condition 25], also used 
below, introduces another length scale. As, the effective 
width of the Stern layer, which also affects the time scales 
for electrochemical relaxation. 

Other important surface properties have also been in- 
cluded in mathematical analyses of AC response. For ex- 
ample, several recent studies of blocking electrodes have 
included the effect of a nonzero equilibrium zeta poten- 
tial (away from the point of zero charge) 65, S^J stIIs ^ . 
building on the work of Delacey and White [83 • A 
greater complication is to include Faradaic processes at 
non-blocking electrodes through boundary conditions of 
the Butler- Volmer type Is^ l57l |. as suggested by 
Levich and Frumkin [91|. This approach has been 
followed in various aiialyses of AC response around base 
states of zero [UlTliEililil and nonzero IHil 
steady Faradaic current. Numerical solutions of the 
time-dependent Nernst-Planck equations have also been 
developed for AC response and more general situa- 
tions |89l l9^ l98l , following the original work of Cohen 
and Cooley |99| . 



C. Colloids and Microfluidic Systems 

Diffuse-charge dynamics occurs not only near elec- 
trodes, but also around colloidal particles and in mi- 
crofiuidic systems, where the coupling with fluid flow 
results in time-dependent, nonlinear electrokinetic phe- 
nomena. This review may be the first to unify some of 
the fairly disjoint literatures on diffuse-charge dynamics 
in these areas with the older literature in electrochem- 
istry discussed above. Compared to the latter, more 
sophisticated mathematical analyses are often done in 
colloidal science and electro-microfluidics, starting from 
the Nernst-Planck equations for ion transport and the 
Navier-Stokes equations for fluid mechanics in two or 
three dimensions. On the other hand, wi th the no table 
exception of the Ukrainian school |3 llOCl llOlj , less 
attention is paid to surface properties, and simple bound- 
ary conditions are usually assumed (constant zeta poten- 
tial and complete blocking of ions) which exclude diffuse- 
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charge dynamics. 

This might explain why the material time scale, r^i, 
is emphasized as the primary one for double-layer re- 
laxation around colloidal particles [23, IH, , although 
the mixed time scale, Tc = {L/Xjj)tjj, has come to be 
recognized as controlling bulk-field screening by elec- 
trodes [^imiiil. This thinking can be traced back 
to the seminal work of Dcbye and Falkenhagen 0, [t^ 
on dielectric dispersion in bulk electrolytes, mentioned 
above. In that context, when a background field Eb is 
applied to an electrolyte, the relevant geometrical length 
is the size of the screening cloud around an ion, L = Ad, 
over which a voltage, Ei,Xd, is efi'ectively applied. The 
relevant RC time for the polarization of the screening 
cloud is then, Tc = Xd^d/D = td- The possible role of 
geometry is masked by the presence of only one relevant 
length scale, Xd- 

For colloidal particles, which are usually much larger 
than the double-layer thickness, the second time scale, 
Ta = /D, for bulk diffusion around a particle of ra- 
dius, a, becomes important, especially in strong fields. 
If there is significant surface conduction or the parti- 
cle is conducting, the "RC" time scale, Tc — Xoa/D, 
can also become important. In general, double-layer re- 
laxation is thus sensitive to the size and shape of the 
particle. Although it is largely unknown (and rarely 
cited) in the West, many effects involving non- uniform 
double-layer polarization around colloidal particles have 
been studied under the name, "non-equilibrium electric 
surface phenome na" [3^ . as recently reviewed by S. S. 
Dukhin |3lll0H. 

The colloidal analog of our model problem involving a 
blocking electrochemical cell is that of an ideally polar- 
izable, metal particle in a suddenly applied background 
electric field. This situation has received much less at- 
tention than the usual case of non-conducting particles 
of fixed surface charge density, but it has an interest- 
ing history. The non-uniform polarization of the double 
layer f or a m etal particle was perhaps first described by 
Levich |102| . usin g Helm holtz' capacitor model. Simonov 
and Shilov |l03l Il04j later considered diffuse charge 
and showed that the metal particle acquires an induced 
dipole moment opposite to the field over the time scale, 
Tc = Xna/D, as bulk conduction transfers charge from 
the part of the double-layer facing away from the field to 
the part facing toward the field. The two hemispheres 
may be viewed a s cap acitors coupled through a continu- 
ous bulk resistor |l04j , as in the RC circuit model of DC 
electrochemical cells described above. The charging pro- 
cess continues until the redistribution of diffuse-charge 
completely eliminates the normal component of the elec- 
tric field, responsible for charging the double layer. 

Diffuse-charge dynamics is important in the context 
of colloids because it affects electrokinetic phenomena. 
In the metal-sphere example, the remaining tangential 
component of the field interacts with the non-uniform 
induced diffuse charge (and zeta potential) to cause non- 
linear electro-osmotic flows |lO0l Il05l | , which causes hy- 



drodynamic interactions between colloidal particles. Al- 
though these fiows have little effect on the electroph oresis 
of ch arged polarizable particles in uniform DC fields |l06l 
Il07| , they significantly affect dielectrophoresis in nonuni- 
form AC fields 108.. .109.] , where the time-dependence of 
double-layer relaxation also plays an important role. 

These developments followed from pioneering stud- 
ies of S. S. Dukhin, B . V. Deryagin, and collabora- 
tors [H, 13 llld llllj on the effects of surface con- 
duction and concentration gradients on electrical polar- 
ization and electro-osmotic flows around highly charged 
non-conducting particles, which was also extended to po- 
larizable particles |lOCl| |. (Similar ideas were also pur- 
sued later in the W e stern liter ature, with some new 
resul ts El Inl ITll ITTl ins.) Earlier stiU, Biker- 
man presented the original theory of sur- 
face conduction in the double layer, and Overbeek |119| 
first calculated in detail the effect of non-equilibrium 
double-layer polarization on electrophoresis. 

Diffuse-charge dynamics has begun to be exploited in 
microfluidic devices, albeit without the benefit of the 
prior literature in electrochemistry and colloidal science 
discussed above. In a series of recent papers, Ramos and 
collaborators have predicted and observed "AC electro- 
osmosis" at a pair of blocking micro-electrodes 0, 0, 
0, 0, 0- Their simple explanation of double-layer dy- 
namics P, , supported by a mathematical analysis of 
AC response in two dimensions 0], is similar to that 
of Si mono v and Shilov for a metal particle in an AC 
field |l04j . and the resulting electro-osmotic flow is of 
the type described by Gamanov et al. for metal parti- 
cles |l05j . An important difference, however, is that AC 
electro-osmosis occurs at fixed micro-electrodes, whose 
potentials are controlled externally, as opposed to free 
colloidal particles. Ajdari has proposed a similar 
means of pumping liquids using AC voltages applied at 
an array of micro-electrodes, where broken symmetries 
in surface geometry or chemistry generally lead to net 
pumping past the array, as observed in subsequent exper- 
iments OS ^1 • These are all examples of the general 
principle of "induced-charge electro-osmosis" 120, 121i|, 
where diffuse-charge dynamics at polarizable surfaces 
(not necessarily electrodes) is used to drive micro-flows 
with AC or DC forcing. Clearly, the full range of possi- 
ble microfluidic applications of time-dependent nonlinear 
electrokinetics has yet to be explored. 



D. The Limit of Thin Double Layers 

All of the analytical studies cited above that go be- 
yond linear response (and most that do not) are based 
on the thin-double-layer approximation, Xd L. In this 
limit, the bulk electrolyte remains quasi-neutral, and the 
double layer remains in thermal quasi-equilibrium, even 
with time dependent forcing (slower than td) Ha, 
IIMUMI^' The same limit also justifles the general no- 
tion of circuit models for the diffuse part of the double 
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layer and, in the absence of concentration gradients, the 
neutral bulk region. 

As first shown by Grafov and Chernenko |l22l 

Il23l | in the Sovi et Union and by Newman |l24j and 
Macgillivray |125| in the United States, the thin double- 
layer approximation for electrochemical cells can be given 
"firm" (but not necessarily "rigorous") mathematical jus- 
tification by the method of matched asymptotic expan- 
sions |33,|40,|41| in the small parameter, e — Xd/L- For 
steady Faradaic conduction, the usual leading-order ap- 
proximation involves a neutral bulk with charged bound- 
ary layers of 0(e) dimensionless width, which has since 
been established rigo r ously in a number of studies by 
mathematicians [iSH EH EH EH EH EHl- The 
standard asymptotic approximation breaks down, how- 
ever, near Nernst's diffusion-limited current, where the 
concentrati on a t the cathode vanishes. At the limit- 
ing current Il32'|, the boundary la yer ex pands to 0{e^/^) 
width, while at still larger currents |l3,3j , a layer of "space 
charge" extends out to 0(1) distances into the bulk re- 
gion, although the effect of realistic boundary condi- 
tions (Faradaic processes, compact layer, etc.) remains 
to be studied in these exotic regimes. Matched asymp- 
totic expansions are also beginning to be used for time- 
dependent electrochemical problems with Faradaic pro- 
cesses IHH^ below the limiting current. 

Perhaps because it originated in fluid mechanics [4ll| . 
the method of matched asymptotic expansions has been 
used extensively in colloidal science and microfluidics 
34, 111, 112, 113, 114, 115,, 134, 13^, albeit with vary- 
ing degrees of mathematical rigor. In any case, the ad- 
vantages of the technique are (i) to justify the assumption 
of equilibrium structure for the double layers (at leading 
order), regardless of transport processes in the neutral 
bulk, and (ii) to view the double layers as infinitely thin 
at the bulk length scale, which is particularly useful in 
multidimensional problems. For statics or dynamics at 
the bulk diffusion time, it is usually possible to construct 
uniformly valid approximations by adding the inner and 
outer solutions and subtracting the overlap. 

The thin double layer approximation is "asymptotic" 
as e — > 0, which means that the ratio of the approxima- 
tion to the exact solution approaches unity for sufficiently 
small e, with all other parameters held fixed. For any 
fixed e > (no matter how small), however, the approx- 
imation breaks down at sufficiently large voltages. The 
general criterion 



^ cosh ( ^ 
a \2kT 



< 1 



(7) 



is often quoted for the validity of Smoluchowski's formula 
for the electrophoretic mobility of a thin-double -laye r 
particle |27j . as justified by numerical calculations |l36j . 
This is related to S. S. Duhkin's seminal work on double - 
layer distortion around a spherical particle [33.lllllll35l| : 
In the case of highly charged particles, C, ^ kT/ ze, the 
"Dukhin number" Du (which he called "Rel") controls 
corrections to the thin-double-layer limit, Du — 0. 



The Dukhin number is defined as the ratio of the 
double-layer surface conductivity, a^, to the bulk con- 
ductivity, (Tf,, in Eq. per geometrical length, a: 
Du = as/(TbO.- Although its effect on electrophoresis was 
first explored in detail by Dukhin, the same dimension- 
less group was defined a few decades earlier by Biker- 
man | 118| . who also realized that it would play a funda- 
mental role in electrokinetic phenomena. In a symmet- 
ric binary electrolyte with equal diffusivities, the Dukhin 
number can be put in the simple form. 



Du = 



2A£,(1 + to) r , f zeC 
4A£i(l + m) . ,2 f 



- 1 



■ sinh 



V 4fcT 



where 



kT 

ze 



2e 



(8) 



(9) 



is a dimensionless number giving the relative importance 
of electro-osmosis compared with electro-migration and 
diffusion in surface conduction, and rj is the viscosity. 
This form is due to Deryagin and D ukhin |l37l , who gen- 
eralized Bikerman's original results |ll(Hll7i | to account 
for electro-osmotic surface conductance (m > 0). For 
Du <^ 1 the double-layer remains in its equilibrium state 
at constant zeta potential, but for Du 3> 1 it becomes 
distorted as surface conduction draws current lines into 
the double layer. For a detailed pedagogical discussion, 
we refer to Lyklema |29l] . 

It is interesting to note that (at least at large zeta 
potentials) the Dukhin number is similar to the ratio of 
the effective RC time, Tc(C), away from the point of zero 
charge {( ^ 0) to the bulk diffusion time, Tq: 



. ^ cosh 



(10) 



where we have used Eqs. ©-(jSJ. Moreover, the usual 
condition for the validity of the thin double-layer 
approximation in quasi-steady electrokinetic problems is 
also a statement about time scales: Tc(C) <C tl- When 
this condition is violated, the usual RC charging dynam- 
ics is slowed down so much by nonlinearity that bulk dif- 
fusion may complicate the picture. Whether this does in 
fact occur depends on if the nonlinearity is strong enough 
to cause significant concentration depletion in the bulk 
for a given geometry and forcing. Understanding this 
issue requires going beyond leading order in asymptotic 
analysis, which is not trivial. 

In spite of extensive work on the asymptotic theory 
of diffuse-charge dynamics, difficult open questions re- 
main. The leading-order thin-double-layer approxima- 
tion is well understood in many cases, but higher-order 
corrections have been calculated in only a few heroic in- 
stances, such as th e asy mptotic analysis of diffusiophore- 
sis by Prieve et al. |ll4j . Moreover, such detailed analysis 
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has mostly (if not exclusively) been done for quasi-steady 
problems. For time-dependent problems of double-layer 
charging, it seems that higher-order terms in uniformly- 
valid matched asymptotic expansions have never been 
calculated. 

Even the leading-order behavior is poorly understood 
when the induced zeta potential is large enough to violate 
the condition Q. In that case, the effective Dukhin num- 
ber varies with time, as the total zeta potential evolves 
in time and space. On the other hand, the Russian lit- 
erature on non-equilibrium electro-surface phenomena at 
large Du mostly pertains to highly charged particles in 
weak fields, where the constant equilibrium zeta potential 
is large, but the time-dependent induced zeta potential 
is small. 

Below, we begin to explore these issues in the much 
simpler context of a one-dimensional problem involving 
parallel-plate electrodes, which excludes surface conduc- 
tion and electro-osmotic flow. We shall see that this re- 
quires extending standard boundary-layer theory, which 
deals with multiple length scales, to account for simul- 
taneous multiple time scales. Before examining the non- 
linear theory, however, we state the mathematical model 
and study its exact solution in the linear limit of small 
potentials. 



III. THE BASIC MATHEMATICAL PROBLEM 

As the simplest problem retaining the essential fea- 
tures of diffuse-charge dynamics, we consider a dilute, 
completely dissociated z : z electrolyte, limited by two 
parallel, planar, blocking electrodes at X = zLL. We 
describe the concentrations of the charged ions by con- 
tinuum fields, C±{X, r), which satisfy the Nernst-Planck 
equations. 



dC± 



d 



(11) 



(without generation/recombination reactions), where 
$ is the electrostatic potential, which describes the 
Coulomb interaction in a mean-field approximation. For 
simplicity we assume that the diffusion coefficients of the 
two ionic species are equal to the same constant, D, 
and obtain the mobility, fi, from the Einstein relation, 
/X = D/kT. The total ionic charge density, pe, controls 
the spatial variation of the potential, through Pois- 
son's equation. 



'dX^ 



Pe = ze(C+ — C_) 



(12) 



where e is the dielectric permittivity of the solvent, as- 
sumed to be a constant. 

As described above, we focus on "ideally polarizable" 
or "completely blocking" electrodes without Faradaic 
processes, so the ionic fluxes have to vanish there: 



r,dC± zeD 9$ 



0, for X = ±L. (13) 



The Faradaic current density, J = ze{F+ — F_), also 
vanishes at the electrodes, although it can be nonzero 
elsewhere as diffuse charge moves around inside the cell. 
We also take into account the intrinsic capacitance of the 
electrode surface through a mixed boundary condition for 
the potential IsR Iflq . The surface capacitance rnay 
represent a Stern layer of polarized solvent molecules W| 
and/or a dielectric coating on the electrode ^3]. If V±{t) 
is the external potential imposed by the external circuit 
on the electrode at X = zLL, then we assume 



$ = ^±TAs7tt7, atX 
oX 



(14) 



where Ag is an effective thickness for the compact part 
of the double layer. For a simple dielectric layer, this 
is equal to its actual thickness times the ratio, s/es, of 
dielectric constants of the solvent and the Stern layer, £5. 

In order to study nonlinear effects and avoid imposing 
a time scale, we consider the response to a step in voltage 
(a suddenly applied DC voltage), rather than the usual 
case of weak AC forcing. For times r < 0, no voltage is 
applied, and we assume no spontaneous charge accumu- 
lation at the electrodes. The initial ionic concentrations 
are uniform, C±{X,t < 0) = Cq. For r > 0, a volt- 
age difference 2V is applied between the two electrodes, 
V±{t > 0) = ±F, and we solve for the evolution of the 
concentrations and the potential. As t ^ 00, the bulk 
electric field at the center, |£'(0,t)| — d^/dX, decays 
from its initial value, V/L, to zero, due to screening by 
diffuse charge which is transferred from the right side of 
the cell (0 < A < 1) to the left (-1 < A < 0). The 
relaxation is complete when the Faradaic current decays 
to zero in steady state, from its initial uniform value, 
J(A,0) = Jo = ~atV/L = -2{zefCoDV/kTL. 



IV. 



LINEAR DYNAMICS 



A. Transform Solution for Arbitrary \d, As, L 

For applied potentials much smaller than the ther- 
mal voltage, V <^ ksT/ze^ the equations can be lin- 
earized, C± = Co + 5C±, so that the ionic charge den- 
sity, Pe = ze{C+ — C-) = ze(6C+ — obeys the 
Debye-Falkenhagen equation T4j i 



D dr dX^ 



K Pe 



(15) 



where k = Xj^ is the inverse screening length. This 
equation can also be written as a conservation law. 



dpe_ _ dJe 
57 ~ ^dX 



(16) 



in terms of the linearized total ionic electrical current, 

(17) 



"^^ ^dX ^^'dX 
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which vanishes at the blocking electrodes, X = ±L. 

To solve the model problem, which involves a step- 
potential in time, it is convenient to use Laplace trans- 
forms, defined by 



fiS) = / dr e-^V(r) 



(18) 



As Pe{X) = for T < 0, the Laplace transforms of 
Eqs. ((121) and are 



where 



s 



"2- — _ I ,,2 

D 



k{S)- = — + K 



(19) 
(20) 

(21) 



The general antisymmetric solution to Eq. (|19|) is, 

PeiX,S) ^ Asmh{kX) (22) 

for some constant A{S), which, substituting into Eq. H2U|) 
and integrating, yields, 

~e^,^{X,S)^^cosh{kX) + B, (23) 

where the constant B{S), determined by Je(±L, S) — 0, 
is given by 



B = Ak cosh{KsL){K ^ - k 



(24) 



Integrating Eq. I|23|l again and enforcing antisymmetry 
yields the Laplace transform of the potential, 

r,^ ,cosh(fcL) /sinh(fcX) kSX \ , , 
$(X, S) = -A ' ( —^Y^ + -y^. ) (25) 



(26) 



efc2 \^cosh(fcL) K^D 
The remaining constant, 

-k^eVS-^^c<MkL) 



kSL f -I I Xs 



tanh(A:L) + Asfc + (1 , ^ 



is determined by the Stern-layer boundary condition, 
Eq. CH). 



B. Long-time Exponential Relaxation 



1. Diffuse Charge Density at a Surface 

Let us focus on one quantity, for example, the Laplace 
transform of the charge density at the anode, pdX — 
L,S). The exact formula is 



Pe{L,S) =ylsinh(fcL), 



(27) 



which is difficult to invert analytically. (Keep in mind 
that k depends on S.) For times much longer than the 
Debye time, we consider the limit, S <^ k^D, in which 
the Laplace transform takes the much simpler asymptotic 
form. 



where 



l+TpS 



1 + kXs coth(KL) 



(28) 



(29) 



and 



L 



coth(KL) (1 



3As ' 
2L , 



1 + K,\s coth(KL) 



(30) 

Since the Laplace transform of l—exp(—T/ro) is S' ^/(H- 
5to), this result clearly shows that the buildup of the 
charged screening layer occurs exponentially over a char- 
acteristic response time given by Eq. 130|l . which is of 
order, L/kD — \dL/D = Tc, for both thin and thick 
double layers. Corrections introduce other mixed scales 
involving the Stern length, such as XsL/D and XsXd/D, 
as well as the Debye time, Xj~i/D. 

Note that the same time scale can also be seen in the 
linear response to a weak oscillatory potential, V± = 
±yRe(e*"'^), which naturally leads to 



(31) 



Pe{L,T) ^ KpRe f Y 



+ iiUTp 



for frequencies well below the Debye frequency, lo <C 
lud ~ D/Xjj. Similar results for AC response near the 
point of zero charge have been obtained by many au- 
thors, as cited above. The characteristic frequency, coc — 
1/tc « D/XdL, also arises the context of AC electro- 
osmotic fluid pumping near micro-electrodes P, 0] , be- 
cause diffuse-layer charging controls the time-dependence 
of the effect. 



There is a great deal of information about transients 
in the Laplace transform of exact solution to the linear 
problem. For times much smaller than the Debye time, 
T <C T£) = Xj^/D (or S ^ K^D), there is no significant 
response, so we are mainly interested in the response at 
longer times, t ^ td (or S <C n^D). There are many 
ways to see that this is generally an exponential relax- 
ation dominated by the mixed time scale discussed above, 
Tc — XdL/D, although several other time scales allowed 
by dimensional analysis also play a role. 



2. Total Diffuse Charge in an Interface 

We now show that the same form of long-time exponen- 
tial relaxation, with a somewhat different characteristic 
time, also holds for other quantities, such as the total 
diffuse charge near the cathode, 

Q{t) = J°^p{X,t)dX, (32) 
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which plays a central role in the nonlinear analysis be- 
low. In the limit of thin double layers, this is simply the 
total interfacial charge (per unit area) of the diffuse part 
of the double layer. Here we consider the total diffuse 
charge near a surface more generally, even when the De- 
bye screening length is much larger than the electrode 
separation. In the latter case, the concept of an "inter- 
face" is not well defined, since the two sides of the cell 
interact very strongly, but we can still study the overall 
separation of diffuse charge caused by the applied volt- 
age. 

Using Eqs. H2UI) and (|23|l . the Laplace transform of the 
total cathodic charge is. 



Q{S) = Afc~i [1 - cosh(fcL)] , 



(33) 



Once again, this is difhcult to invert analytically, so we 
focus on the long-time limit. 



QiS) 



KqS-^ 
1 + tq5' 



for S ^ K^Z?, where 



Kr 



eVk [1 - sech(fcL)] 
tanh(KL) -I- kXs 



(34) 



(35) 



and 




tanh(KL) -|- kAs 
sech(Ki) tanh(KL) 



2 [1 - sech(KL)] 



In the limit of thin double layers, the same basic time 
scale, Tc = L/kD — \dL/D, arises as in the case of 
the surface charge density. A subtle observation is that 
the relaxation of the total interfacial charge, although 
still exponential, has a somewhat different time scale as 
a function of e = \d/L and 5 = Xs/\d- (See Fig. [21 
below.) This apparently new result shows that charging 
dynamics has a nontrivial dependence on time and space, 
even for very weak potentials. 



V. DIMENSIONLESS FORMULATION AND 
NUMERICAL SOLUTION 

A. Basic Equations 

In preparation for analysis of the full, nonlinear prob- 
lem, we cast it in a dimensionless form using L as the ref- 
erence length scale and Tc = XdL/D as the reference time 
scale, as motivated by the linear theory. Time and space 
are then represented hy t — tD/X^L and x — X/L^ 
and the problem is reformulated through reduced vari- 
ables: c = (C+ -f C^)I2Cq for the local salt concentra- 
tion, p = (C+ — C_)/2Co = Pe/{'2.Coze) for the charge 
density, and = ze^/ksT for the electrostatic potential. 




1.0 0.5 1.0 1.5 2 

e 



0.5 1.0 1.5 2.( 

e 



FIG. 2; Analytical results for the exponential relaxation 
time from the linear theory for weak applied potentials, 
{V <^ kT/ze). The time scale for relaxation of the surface 
diffuse-charge density, tp, from Eq. 1451 is shown in (a) and 
that of the total interfacial (half-cell) diffuse charge, tq , from 
Eq. 15UII in (b). In each case, the charging time, scaled to 
Tc = L/kD = \dL/D, is plotted versus the dimensionless 
diffuse-layer thickness, e = \d/L, for different dimensionless 
Stern-layer thicknesses, 5 — As/Au ~ 0,0.1, 1, 10 (solid, dot, 
dash, and dot-dash lines, respectively). 



The solution is determined by only three dimensionless 
parameters: v = zeV/ksT, the ratio of the applied volt- 
age to the thermal voltage, e = Xd/L, the ratio of the 
Debye length to the system size, and 6 = Ag/An, the 
ratio of the Stern length to the Debye length [93 . 

With these definitions, the dimensionless equations for 
— 1 < a; < 1 and t > are 



(37) 
(38) 
(39) 

(40) 
(41) 
(42) 



dc 


d 


( dc 


dcf) 


dt 


^ 'Tx 


\dx ~ 


'"Tx 


dp 




(dp 


d(t) 


at 


^ 'Tx 


\dx 


dx 










dx"^ 


= p 







with boundary conditions at x = ±1, 



dc 
dx 

dp 

dx 
V — Se 



dx 



' dx 



dx 








and initial conditions, c{x,0) — 1, p{x,0) — 0, and 
(f){x, 0) = V X. Note that the limit of a negligible screening 
length, e ^ 0, is singular because it is impossible to sat- 
isfy all the boundary conditions when e = 0. Physically, 
this corresponds to the limit of exact charge neutrality, 
p = 0, which is always violated to some degree at elec- 
trochemical interfaces. 
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The total diffuse charge near the cathode is 



= y p{x,t)dx, 



(43) 



scaled to 2zeCoi. The dimensionless Faradaic current 
density is, 



dp d(j) 
ox ox 



(44) 



scaled to 2zeCQD / L (Nernst's diffusion-limited cur- 
rent 113). 



B. Time Scales for Linear Response 

The time scale for exponential relaxation of the surface 
charge density in the linear theory above, Eq. 130|) . has 
the dimensionless form, 



(1 + 3Se/2) coth(e-i) - (5csch^(e-i)/2 - e 
1 -|-(5 coth(e-i) 



(45) 



As shown in Fig. [21 a), this formula shows that for a wide 
range of diffuse and Stern layer thicknesses, the basic 
time scale is always roughly of order, X^L/D, since tc is 
of order 1. In the limit of a thin diffuse double layer, the 
dimensionless time scale has the form. 



3(5-2 



1 + 5 \2{1 + S) 



e + 0{e~' ). (46) 



with exponentially small errors. In the limit of a thin 
Stern layer, the time scale becomes 

tp = coth(e"^) - e + 



[5ecoth(e"^) - 2coth(e"^) - csch^(e 
+ 0{S^). 



(47) 



For simultaneously thin Stern and diffuse layers, we ob- 
tain the simple result. 



tr. 



1-e-S 



(48) 



which, as in Fig. [2Ia), shows that increasing either e = 
Xd/L or 6 = Xs/Xd tends to reduce the charging time in 
this limit, compared to the leading-order value, XdL/D. 
Putting the units back, this expression can be written as. 



XpL 
D 



X2 

D 



D 



(49) 



for As ^ Ad ^ L, which clearly shows the Debye time, 
X^/D, appearing only as a small perturbation of the in- 
termediate time scale, X]jL/D, for the relaxation of the 
ceh. 



Similar results hold for the relaxation time for the total 
half-cell charge, Eq. (|36|) . which has the dimensionless 
form, 



l-f isech2(e-i) + 
tanh(e-i) + S 



€ sech(e ^)tanh(e ^) 

2 " 



2 [1 - sech(e-i) 



(50) 

For thin double layers, we obtain the same leading-order 
behavior, 



1 



1 + S 



1 - 26 



2{l + 5) 



+ 0{e~^- ) 



(51) 



although the correction term is somewhat different for 
thick diffuse layers. For simultaneously thin diffuse and 
Stern layers, the dimensionless relaxation time for the 
total charge becomes. 



1 



(52) 



For a detailed summary of how the two time scales, tp 
and tg, depend on the parameters, e and (5, see Figure El 
(a) and (b), respectively. 



C. Numerical Solution 

Our dimensionless model problem, stated in Section 

IV Al is straightforward to solve numerically using finite 
differences, at least if e is not too small. (Ironically, as 
shown below, analytical progress is much easier in this 
singular limit.) To resolve the boundary layer where 
the gradient is large, a variable size mesh is used, along 
with second-order-accurate differencing that accounts for 
the variable grid sizes. The third-order Adams-Bashforth 
method is used in time. The number of the grid points 
and the ratio of the smallest to largest grid size are var- 
ied depending on the values of e and v. The numerical 
convergence is verified though multiple runs of different 
resolutions, and as a result, up to 1024 points are used 
in calculations for higher v. 

To maximize the importance of diffuse charge, we first 
consider a rather larger value of e, even for a micro- 
electrochemical system, e = 0.05, say for Ad = 5 nm and 
L — 0.1/im. The Stern length is always of molecular di- 
mensions, so we choose As = 5 A, and thus 6 = 0.1. The 
time evolutions of the charge and potential are shown in 
Fig. Ofor V — 0.1 and v — 1. At room temperature (and 
2 = 1), these voltages correspond to V — 2.5 mV and 

V — 25 mV, respectively, which, when transferred to the 
diffuse layer after screening give maximum electric fields 
of order 10 V/pm. 

The current, j, and the total cathodic diffuse charge, q, 
are plotted versus time, t, in Fig. 0] for apphed voltages, 
V = 1,2,3, and 4. In all cases, the linearization is accu- 
rate at early times (t < 1) since the dimensionless voltage 
across the diffuse layer remains small (< 1). For v = 1, 
the linear approximation is reasonable for all times, but 
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for somewhat larger voltages, v = 2, 3, and 4, the rel- 
ative error becomes unacceptable at long times, i > 1. 
Not only is the limiting value of the total charge signifi- 
cantly underestimated, but the dynamics also continues 
for a longer time, with a qualitatively different charging 
profile. The largest applied voltage, v — A, shows this 
effect most clearly, as there is a secondary relaxation at 
a much larger time scale of order t — 1/e = 20. Un- 
like the other cases, which display the expected steady 
increase in charge of an RC circuit, for w > 4 the total 
charge quickly reaches a maximum value, after the initial 
RC charging process, and then slowly decays toward its 
limiting value. 

We are not aware of any previous theoretical prediction 
of such a non-monotonic charging profile, so it is a major 
focus of this work (in sections IVlll and lVlll() . It is rem- 
iniscent of the Warburg impedance due to bulk diffusion 
of current-carrying ions at the time scale, tl = L-/D, 
or 1/e in our units, in (linear) response to Faradaic pro- 
cesses, which consume or produce them at an electrode. 
Here, however, there are no Faradaic processes, so any 
such bulk diffusion must be related to the adsorption 
or desorption of ions in the diffuse part of the double 
layer. Moreover, the over-relaxation of the charge den- 
sity is part of the non-linear response to a large applied 
voltage, so it will require more sophisticated analytical 



methods. 



VI. WEAKLY NONLINEAR DYNAMICS 
A. Asymptotic Analysis for Thin Double Layers 

The remarkable robustness of the charging time well 
into the nonlinear regime (at least for the primary relax- 
ation phase) can be predicted analytically using matched 
asymptotic expansions in the singular limit of thin dou- 
ble layers, e = Xd/L <C 1. Most (if not all) previ- 
ous studies of time-dependent problems using asymp- 
totic analysis have scaled time to the diffusion time, 
tl — L^/D. In this section, we will see how the cor- 
rect charging time scale, Tc = XdL/D, arises systemati- 
cally from asymptotic matching at leading order. We also 
consider, apparently for the first time, the general case of 
arbitrary voltage, v — zeV/kBT, and Stern-layer thick- 
ness, S = As/ Ad, with a time-dependent zeta-potential 
(i.e. potential drop over the diffuse layer). We also study 
higher-order corrections, which involve some bulk diffu- 
sion at the time scale, r^. 

As usual, matched asymptotic expansions only pro- 
duce a series of "asymptotic" approximations to the so- 
lution, in the sense that higher terms in the expansions 



12 




1.2 

0.8 

0.4 
0.0 

I 

1.4 
1.2 

1.0 
0.8 



' 1 ' 1 — 




1 1 1 ■ 1 


(b) ■ 







2 3 
t 



(C) 







10 15 20 25 30 
t 



FIG. 4: (a) The dimensionless current density, j(t) (in units 
of 2zeCoD / L), and (b) the dimensionless total diffuse charge 
on the cathodic side of the cell, q{t) (in units of 2zeCoL), 
scaled to qo = v/{l + 5), versus dimensionless time, t (in units 
of Tc = \dL/D). Numerical results for dimensionless volt- 
ages, II = 1 (dot), 2 (dash), 3 (dot-dash), and 4 (dot-dot-dot- 
dash) are compared with linear dynamics in the thin double- 
layer limit: q{t)/qa ~ 1-6-'^+*'' and j(t)/v ~ e-(^+*)* (solid 
lines) as v,e 0. The breakdown of linear theory for v > 1 
is highlighted in (c), where the data in (b) is replotted for 
longer times. 



vanish more quickly than the leading terms as e ^ 0, 
with the other parameters, v and 6, held fixed at ar- 
bitrary values. For any fixed e > (no matter how 
small), there could be e-dependent restrictions on v and 
S for various truncated expansions to produce accurate 
approximations. We refer to the regime where such con- 
ditions hold as "weakly nonlinear" , as opposed to the 
"strongly nonlinear" regime where the asymptotic expan- 



sions break down (described below in section IVIIip . 

B. Outer and Inner Expansions 

We begin by seeking regular asymptotic expansions 
(denoted by a bar accent) in the bulk "outer" region, 
e.g. 



c{x, t) ~ c{x, f) = Co + e ci + C2 + 



(53) 



Substituting such expansions into Eqs. (|T7|l - l|S^ and 
equating terms order by order yields a hierarchy of 
partial-differential equations. At leading order in e, we 
find that the bulk concentration does not vary in time, 
Co = 1, simply because the charging time scale, Tc, is 
much smaller than the bulk diffusion time scale, t^. The 
leading-order potential is linear. 



00 = jo{t)x, 



(54) 



where jo(0) = Since the leading-order bulk concentra- 
tion is uniform, jo{t) is the leading-order current density. 
The leading-order charge density, 



P2 



(55) 



vanishes because the leading-order potential, Eq. (|54|l . is 
harmonic, although at next order, 0{e'\ a nonzero bulk 
charge density, pa, arises due to concentration polariza- 
tion. (See below.) These arguments justify the usual as- 
sumption of bulk electroneutrality to high accuracy, even 
during intcrfacial charging, as long as the dynamics are 
"weakly nonlinear" . 

The regular outer approximations must be matched 
with singular "inner" approximations in the boundary 
layers. The problem has the following symmetries about 
the origin, 



c(— a;, t) ~ c{x, t) 
p{~x,t) = -p{x,t) 
cj){-x,t) = -cj){x,t) 



(56) 



so we consider only the boundary layer at the cathode, 
X — —1, by transforming the equations to the inner co- 
ordinate, y — [x + 



(57) 



dc 




f d~c 


.d4> 




dy 




dy 


dp 


= A 


( dp 






dy 


\dy' 


dy 












= ~P 







(58) 



(59) 



This scaling removes the singular perturbation in Pois- 
son's equation, so we can seek regular asymptotic ex- 
pansions for the inner approximations (denoted by tilde 
accents), e.g. 



c{x, t) ^ c{y, t) = Co + e ci + C2 + 



(60) 
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Matching with the bulk approximations in space involves 
the usual van Dyke conditions, e.g. 



lim c{y,t) ^ lim c{x,t), 



(61) 



which implies co(oo,i) = co(— l,t), ci(cx),i) = ci(— l,t), 
etc., but we will also have to make sure that the expan- 
sions are properly synchronized in time. In particular, 
we will have to worry about the appearance of multiple 
time scales at different orders. 

Substituting the inner expansions into the rescaled 
equations (|57l) - (|5^ causes the time-dependent terms 
to drop out at leading order. Physically, this quasi- 
equilibrium occurs because the charging time, Tc, is much 
larger than the Debye time, td, characteristic of local dy- 
namics in the boundary layer (at the scale of the Debye 
length, Xd)- As a result, we systematically arrive at clas- 
sical Gouy-Chapman profiles for the equilibrium diffuse 
layer at leading order, 



c± ~ e^"^, Co = cosh ipQ, pq = — sinh-^o 
where the excess voltage relative to the bulk, 
i>iy, t) = 4>{y, t) ~ ^(™1, t) ~ + e V^i + • • • 



(62) 



(63) 



satisfies the Poisson-Boltzmann equation at leading or- 
der. 



dip- 



sinh ■(/'o. 



(64) 



Note that matching implies "00(00; i) = V'i(oOji) — . . 
0. The dimensionless zeta potential, C(t) = 0(0, i), varies 
as the diffuse layer charges. 

After the first integration we apply matching to the 
electric field, 



-(oo,i)^.-(-l,t) 



— — c»,t) = 0, 
dy 



to obtain 



^ = -2 sinh(^o/2). 
dy 

After the second integration, 

^o{y,t) = -4 tanh""^(e 
we are left with a constant, 

i^(i) =logcoth(-Co(i)/4) 



(65) 



(66) 



(67) 



(68) 



to be determined from Co (t) (below) by the Stern bound- 
ary condition at the cathode surface, y = 0, and the 
coupling to the bulk region. The offset parameter, K{t), 
which also appears in the concentration and charge den- 
sity, 

co(y,t) = l + 2csch^{y + K{t)) (69) 
Poiy,t) = 2csch{y + K{t))coth{y + K{t)) (70) 

is quite sensitive to Faradaic reactions [o^, but here we 
focus only on the effect of compact-layer capacitance. 



C. Time-dependent Matching 

It seems we have reached a paradox: Both the bulk and 
the boundary layers are in quasi-equilibrium at leading 
order, and yet there must be some dynamics, if we have 
chosen the proper time scale. The resolution lies in taking 
a closer look at asymptotic matching. Physically, we are 
motivated to consider the dynamics of the total diffuse 
charge, which has the scaling, q(t) ~ eq{t), where 



P(y, t)dy qo + eqi+ q2 + 



(71) 



Taking a time derivative using Eq. (|58|l and applying the 
no-flux boundary condition H41|) . we find 



dq ^. 1 f dp ^ _d(f) 
dt y~*oo e \dy dy 



where we have applied matching to the derivatives (flux 
densities). Substituting the inner and outer expansions 
yields a hierarchy of matching conditions. At leading 
order, we have 



f (t)=3o(t), 



(73) 



which shows that we have chosen the right time scale 
because this is a balance of 0(1) quantities. More- 
over, it can be shown that any other choice of scaling 
would lead to a breakdown of asymptotic matching in 
the limit e ^ 0. (For example, in the analogous Equa- 
tions (42)-(43) of Ref. ilj for small AC potentials, the 
time-dependent term vanishes in this limit, showing that 
the proper scaling was not used.) Therefore, the cor- 
rect charging time scale, Eq. |(2J), in the weakly nonlin- 
ear regime follows systematically from time-dependent 
asymptotic matching at leading order. 

The physical interpretation of Equation H73|l is clear: 
At leading order, the boundary layer acts like a capacitor, 
whose total charge (per unit area) , q, changes in response 
to the transient Faradaic current density, j{t), from the 
bulk. The matching condition can also be understood 
physically as a statement of current continuity across the 
diffuse layer. Substituting Poisson's equation into 
Eq. (|71|l . integrating, and matching the electric field us- 
ing Eq. linnj), we see that the left hand side of Eq. iTTl^ 
is simply the leading-order (dimensionless) displacement 
current density [osl lo^ IqtL lo^ at the cathode surface. 



dqo d d(j)o ■-. 



(74) 



so the matching condition simply reads, jo{t) — ia{t)- 
This transient displacement current exists in the external 
circuit, even if there is no Faradaic current. 
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FIG. 5: Sketch of the equivalent RC circuit for the leading- 
order weakly nonlinear approximation: compact-layer and 
diffuse-layer capacitors in series with a bulk resistor. Al- 
though remarkably robust, the circuit approximation is vio- 
lated by higher-order corrections, especially at large voltages. 



D. Leading-order Dynamics 

Using Eqs. ^ and the integral in Eq. lf7T|) 
can be performed at leading order to obtain the Chap- 
man's formula for the total diffuse charge, 

qo = -2 sinh(Co/2). (75) 

The Stern boundary condition, Eq. 142(1 . then yields, 



Co + 25 sinh(Co/2) - jo{t) -v = ^o, 



(76) 



where ^{t) = -v - (j>{-l,t) §o + f^i + ... is the total 
voltage across the compact and diffuse layers. Substi- 
tuting into the matching condition, Eq. 1(7^^(1 . we obtain 
an ordinary, initial- value problem, either for the leading- 
order double-layer voltage 



-Co(*o)^ = *o + ^^ *o(0) = 0, 
or for the leading-order current density 
djo 



Co (jo - w) 



dt 



-Jo, jo(0) = w, 



(77) 



(78) 



where Co(^'o) = dqo/d'^o is the differential capacitance 
for the double layer as a function of its total voltage, 
relative to the potential of zero charge. 

The effective double-layer capacitance is given by 



Co = 



1 



sech(Co/2) + (5 



(79) 



where Co is related to by Eq. (|76|l . A similar formula 
arises in the classical circuit model of Macdonald '63*1 . In- 
deed, the leading-order charging dynamics from asymp- 
totic analysis corresponds exactly to the nonlinear RC 
circuit shown in Fig. We expect, however, that the ad 
hoc circuit approximation cannot describe higher-order 
asymptotic approximations, where the finite thickness of 
the double layer becomes important. 



Linearizing for small voltages, Cq ~ + 5), we ob- 
tain the same results as before in the limit e ^ 0, now 
by a completely different method. 



hit) 

qo{t) 



V e 
v{l 



(1+S)t 



-il+S}t 



1 



(80) 
(81) 



As shown in Fig.^for 6 = 0.1, the linearization describes 
the charging dynamics fairly accurately, even for some- 
what large voltages (w w 1), as long as S is not too small. 
One way to understand this is that the total differential 
capacitance satisfies the uniform bounds, 

7^7 = ^o(O) < Coi^o) < Coi^) = (82) 

1 + 

in the linear and nonlinear regimes. Moreover, the lin- 
earization is always accurate at early times (up to i « 1 
or T « Tc) for any applied voltage, as long as the initial 
zeta potential (or diffuse charge) is small. This is also 
clearly seen in Fig. ^ 

The dynamical equation H77|) or (|78|l is first-order and 
separable, so its exact solution is easily expressed in in- 
tegral form, 



^oit)=joit)~v = -F-\t) 



where 



F(z) 



Co (u) du 



(83) 



(84) 



The integral can be evaluated numerically and the to- 
tal charge recovered from the Eqs. H75|l and ((76|l . The 
results in Fig. IHl show that the leading-order dynamics 
compares fairly well with the numerical solution to the 
full nonlinear problem for e — 0.05 and 5 = 0.1, at 
least for the decay of the current density, especially at 
early times (t ~ 1). The limiting value of the total dif- 
fuse charge is also approximated much better than in 
the linear theory (Fig. due to the nonlinear differ- 
ential capacitance, Eq. (|79|l . For large voltages {v > 1), 
however, total charge shows some secondary dynamics at 
longer time scales {t ^ 1), which is not fully captured 
by the leading-order asymptotic approximation (or the 
corresponding circuit model). As we shall see below, this 
can only be understood by considering higher-order terms 
which violate the circuit approximation. 

For moderately large voltages (w « 1), we can expand 
around u = v in the integrand of Eq. H84|) and obtain a 
long-time exponential decay. 



jo(t) = « + *o(i)oce-*/^«(") 



(85) 



as t —f oo. This reveals a (dimensionless) characteristic 
time, tc = Co{v), which is larger than that of the linear 
regime, tc = Cq{0) = 1/(1 + (5), by at most a factor of 
1/(5 {— 11 in our numerical examples). Although this 



1 



factor is non- negligible, the characteristic time, t^, is still 
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FIG. 6: The comparison of the full numerical solution with 
the leading-order asymptotic results: (a) j/u, (b) q/qun (early 
evolution), and (c) q/qnn (long time evolution). The full nu- 
merical solution are shown with dot {v — 1), dot-dash (with 
open squares in (c)) {v = 2) and long dash with open trian- 
gles {v = 3). The leading-order asymptotic results are plotted 
with dash {v = 1), dot-dot-dot-dash (with filled squares in 
(c)) {v = 2) and long dash with filled triangles (« = 3). The 
curves for v — 3 are omitted in (a) and (b) for clarity. The 
solid lines show the linear dynamics in the thin double-layer 
limit. 



the basic time scale, rather than and tu which differ 
from Tc by factors of e, i.e. usually two or more orders of 
magnitude. As the voltage is increased, however, nonlin- 
earity always becomes important, and one of its generic 
effects is to slow down the relaxation process. 

In order to simplify the response function, F{z), and 
other quantities, it is useful to consider the regular limit 
of thin Stern layers, (5 — *■ (taken after the singular limit 
of thin diffuse layers, e — > 0). In this common physi- 
cal regime, where As «C Ad <C i, the following asymp- 
totic expansions can be derived by iteration |40j from 



Eqs. (ESI), (EH, and 

Co *o - 25 sinh + (5^ sinh 



90 



5*0 

-2 sinh + (5sinh*o 
*o 



-(5^ I sinh *o cosh — + sinh^ — 



*o 



(86) 



(87) 



*o 



Co ~ cosh — S cosh '^q 

+S'^ I cosh ^0 cosh + — sinh sinh ^ 



3 • u2 *0 , ^-0 

-— smh — cosh — 
2 2 2 



The response function can then be expanded in somewhat 
simpler (but still nontrivial) integrals. 



Fiz) 



Jo u + V 
in the limit 6—^0. 



cosh(u/2) d-u ^ f^cosh{u)dx 



u + V 



(89) 



E. Uniformly Valid Approximations 

Asymptotic analysis tells us not only the behavior of 
integrated quantities like total charge and voltage, but 
also the complete spatio-temporal profiles of the charge 
density and potential. As usual, uniformly valid approx- 
imations (in space) are constructed by adding the outer 
and inner approximations and subtracting the overlaps. 
Taking advantage of the symmetries in Eq. H56|l , we ob- 
tain the following leading-order approximations: 



0(a;,i) 
c(x, t) 
p{x,t) 



Col ,t]+co[ ,< -1 (91) 



Po 



e 

1 + x 



^ I - Po 1 - — -,t 



(92) 



where the boundary-layer contributions are given by 
Eqs. ^7^-^^ and Eq. (ESJ, which express the effect of 
the compact layer. The time dependence of the leading- 
order approximations is entirely determined by the bulk 
current density, jo{t), or the double-layer voltage, ^o{t), 
via Eqs. and 

As shown in Fig. [7| the time-dependent approxima- 
tions for (/) and p are in excellent agreement with our 
numerical results well into the nonlinear regime {v = 1), 
even for a fairly large boundary-layer thickness, e — 0.05. 
The charge density clearly shows the expected separa- 
tion into three regions: a neutral bulk with two charged 
boundary layers of 0(e) width. On the other hand, for 
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FIG. 7: The potential (a), charge density (b) and concentration (c) at t = 1 for a large dimensionless voltage, v = 1, 
with e = 0.05 and S = 0.1. The full numerical solution (dashed lines) is compared with the leading-order uniformly valid 
approximation (dotted lines), Eqs. I|9()|l - H92|l . The analytical approximations are almost indistinguishable from the numerical 
solutions for (/> and p, but not for c, which shows errors of a few percent (< e) just outside the double layers. 



the same parameters, the leading-order approximation of 
c is not nearly as good. As expected, the concentration 
exhibits a homogeneous bulk region and two inhomoge- 
neous boundary layers of 0(e) width, which are fairly 
well described, but there are also intermediate regions of 
depleted concentration extending far into the bulk, which 
are not captured at leading order. 



VII. HIGHER-ORDER EFFECTS 

A. Neutral-Salt Adsorption by the Double Layer 

We have seen that each diffuse-charge layer acquires an 
excess salt concentration relative to the outer region. At 
leading order, however, there is no sign of how the extra 
ions got there. This paradox, which also applies to cir- 
cuit models, is apparent from symmetry alone, Eq. H56|l 
— Diffuse charge near the cathode grows by bulk elec- 



tromigration, which creates equal and opposite diffuse 
charge near the anode. In contrast, the excess concen- 
tration is the same in both double layers, so it can only 
arrive there by diffusion of neutral electrolyte from the 
bulk, which is excluded at leading order. 

The key to understanding higher-order terms, there- 
fore, is the total excess concentration per unit surface 
area in (say) the cathodic diffuse layer, w{t) — €w{t), 
where 

/•oo 

w{t)^ / [£{y,t)-co{-l,t)]dy = wo{t)+£Wl{t) + ■■■ 
(93) 

is analogous to the scaled total diffuse charge, q{t). (Note 
that co{—l,t) = 1 in our model problem, but Equation 
(|93ll is more general.) We proceed with matching in the 
same manner as above. Taking a time derivative using 
Eq. H57|) and applying the no-flux boundary condition 
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(|^ . we find 

dw Y ^ \ '^'^ Y ( -^^ 

dt y^oo e \dy ^ dy j x^-i \dx ^ dx 

(94) 

Substituting the inner and outer expansions yields an- 
other hierarchy of matching conditions. At leading order, 
we have 

dwQ 1 dwf) dci 

which, unlike Eq. (|73|l . involves a new time variable, 

t = et^ — ^—, (96) 

scaled to the bulk diffusion time, tl — L? /D. Physically, 
this matching condition simply expresses mass conserva- 
tion: The (zeroth order) excess concentration in the dif- 
fuse layer varies in response to the (first order) diffusive 
flux from the bulk. 

In Eq. (|93|l . the left-hand side is given by the leading- 
order inner approximation calculated above. Substitut- 
ing Eq. (innj into Eq. JHSJ, integrating, and using Eq. 
yields, 



dm , . 
dt ^ ' 



2— coth K{t) 
dt ^ ' 

2-^cosh%) 
dt 2 

go(^) <K^{t) 
2 dt 



(97) 



where we have used the identity, cosh2z = 
— coth log tanhz. Recall that the leading-order zeta po- 
tential, Co(i), is related via Eq. (|76|l to the leading-order 
bulk current density, jo{t), or interfacial voltage, 4*0 (i), 
given by Eqs. ||S2Jl and (jH^. Integrating Eq. and 
requiring wq ~ for (^q — 0, we also obtain a simple 
expression for the excess concentration. 



the "Donnan effect" with refere nce to Donnan' s general 
papers on membrane equilibria jl38l Il39l Il4(| and de- 
scribes how it is used to infer surface areas from exper- 
imental measurements of concentration variations upon 
charging ( "salt sieving" ) j2^ .141.] . In the present case of 
diffuse-layer adsorption, the following argument is given: 
Since the equilibrium co-ion concentration in the double 
layer near a charged surface is reduced compared to the 
bulk, there must be an excess of co-ions, and hence an 
excess of neutral concentration (and counter-ions) in the 
nearby bulk. 

How can this belief be reconciled with our analytical 
and numerical results, which clearly demonstrate the op- 
posite effect in a model problem? The difference is that 
we consider a finite system with global ion conservation, 
while Lyklema considers an open system into which ions 
are apparently injected. We also explicitly calculate the 
time-dependent response to interfacial charging, while he 
describes the steady state after new ions somehow arrive 
"from infinity" . 

Our conclusion is therefore the opposite: Since the 
equilibrium counter-ion concentration is enhanced in the 
double layer, there must be a depletion of counter-ions, 
and hence a reduction in neutral concentration (and co- 
ions) in the nearby bulk. Gouy's nonlinear theory shows 
that at large voltages the excess of counter ions exceeds 
the reduction in co-ions in the diffuse layer, so this result 
is quite intuitive. 

In real systems, it may be that compact-layer effects, 
such as the specific adsorption of ions on the surface, can 
lead to overall negative adsorption by the double layer. 
According to the Nernst-Planck-Poisson equations, how- 
ever, the average concentration of all ions (regardless of 
species) is always increased in the diffuse part of the 
double layer relative to the bulk in any finite system. 
The associated local depletion of neutral salt has impor- 
tant implications for time-dependent electrokinetic phe- 
nomena at polarizable surfaces (section III Cll since bulk 
concentration gradients alter electric fields and produce 
diffusio-osmotic slip. 



Wq = A sinh^ 



(98) 



which also holds for the static Gouy-Chapman solution. 
Of course, this is another sign that (at leading order) a 
thin double layer stays in quasi-equilibrium, even while 
charging. 



B. The Sign of the Donnan Effect 

Before proceeding to calculate the bulk dynamics, we 
comment on the sign of the excess concentration in the 
diffuse part of the double layer, which corresponds to 
a positive adsorption of neutral salt. In contrast, it is 
commonly believed that double-layer salt adsorption is 
always negative, resulting in an excess neutral concen- 
tration in the nearby bulk electrolyte. Lyklema calls this 



C. Bulk Diffusion at Two Time Scales 

We now proceed to calculate how the bulk concentra- 
tion is depleted in time and space during double-layer 
charging in our model problem. The matching condition, 
Eq. H95|l . seems to contradict the analysis above, since it 
introduces a new time variable, i. However, this is the 
same time scale for the first-order (diffusive) dynamics in 
the bulk, 



dci 



1 dci 
7 'dt 



(99) 



We must solve this equation starting from ci(a;,0) = 
with a time-dependent prescribed flux at a; = — 1 given 
by Eqs. (|95|) and (|97|) . We also enforce symmetry about 
the origin, Eq. (|56|l . 
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The Laplace transform of the solution is: 



ci(x, s) 



/s cosh(x-\/s) I 



sinh(Vs) JO 
e~^*ci{x, t)dt 



e^'*wo{t/e)dt 



(100) 



where wo{t) is determined by Co(i) from Eq. H98(l . The 
prefactor, 



G(.) = 



cosh(a;-ys) 
sinh(ys) ' 



(101) 



is the Laplace transform of G{t), the Green function for 
the diffusion equation, Eq. H99|l . for a sudden unit flux 
of ions at time t = 0^ injected at the boundary: 

Gix,0)^0, ^i~l,t)^S+{i). (102) 

The same Green function also arises in the equivalent 
problem of an initial unit source adjacent to a reflecting 
wall, 

BC 

G{x,0) ^S{x + 1+), _(-l,t) = 0. (103) 
ox 

In this form, the Green function can be obtained by in- 
spection, 



oo 



(104) 



using the method of images. 

Since ci{x, s) is expressed as a product of two Laplace 
transforms, Eq. (|10U|I . the inverse is equal to the convo- 
lution of the two original functions: 



c^{x, t) = - f dp G{x, P -t)^^^' 



(105) 



This form clearly demonstrates that the boundary forcing 
occurs over the fast, charging time, t = t/e, while the 
response described by the Green-function kernel occurs 
over the slow, diffusion time, t. The separation of time 
scales is apparent in the equivalent expression, 

c^{x,t)^~ j\t'G{x,e{t' ~t)) tanh|^^^ Ut') 

(106) 

which can be derived from Eq. H105|) using Eq. H97|) . This 
form shows explicitly how solution for the current at lead- 
ing order, Eq. H83|l . fully determines the bulk concentra- 
tion at first order. 

Before further analysis of the exact solution for ci(x, t), 
we describe its physical significance. The bulk concentra- 
tion at first-order exhibits diffusive relaxation at two dif- 
ferent time scales, t = 0(1) and i = 0(1), or with units, 
r = 0{XdL/D) and r = 0{L^/D), respectively. For 




0.90 



FIG. 8: Weakly nonlinear dynamics for i; — 1, e = 0.05, 5 = 
0.1, showing the effect of bulk diffusion. The concentration 
from the full numerical solution is shown in the half cell (a) 
and in the diffuse layer (b) for t = 0.5 (solid), t = 1 (dot), 
2 (dash), 4 (dot-dash), 8 (dot-dot-dot-dash), and 20 (long 
dash) . 



t = 0(1) and t — 0(e), the initial double-layer charging 
process proceeds without any significant changes in con- 
centration at the bulk length scale, x = 0(1). During this 
phase, each diffuse-charge layer acquires an 0(e) amount 
of excess concentration, given by Eq. H109() . This excess 
concentration has been acquired by a diffusive process, 
which at this time scale corresponds to a bulk diffusion 
layer of 0{y/e) width near each electrode. This implies 
an 0(-\/e) depletion of the neutral salt concentration in 
the bulk diffusion layers. These scaling arguments are 
confirmed by Fig.[7|^c) for u = 1 and e = 0.05, where at 
time t = 1 (or t = e) the diffusion layers are roughly of 
width V^i = \/2e « 0.3. The formation and spreading of 
the diffusion layers is also shown in more detail in Fig.|Sl 



D. Evolution of the Diffusion Layers 

In the previous section, we derived the time-dependent 
outer approximation, Eq. H53|l . to first order. 



c(a;,<) ^ 1 -f- eci(x,t). 



(107) 



which displays dynamics at both the RC time and the 
bulk diffusion time. The result, Eqs. H104|l - H106|) . is fairly 
complicated, so in this section we try to gain some simple 
analytical insight. In the limit e — s- 0, the initial charg- 
ing process at the time scale t — 0(e) is instantaneous, 
and we are left with only the slow relaxation of the bulk 
diffusion layers. Explicitly taking this limit in Eq. 1100|) 
with t = et fixed. 



lim ci{x,t) 



-wo{co)G{x,t), 



(108) 



we see that the slow-scale evolution of the diffusion layers 
is given by the Green function, G{x,t), with a source 



19 



1.000 
0.995 

o 0.990 

0.985 



0.980 







1 


1 


1 








/ /' 


\ 












\ 












\ 










/ / 


\ 






















/ / 












1 / 












1 / 




\ \ 








/ /' 




\,\ 












\\ 






















1 







-1.0 



-0.5 



0.0 

X 



0.5 



1.0 



FIG. 9: Simple approximations of the bulk diffusion layers 
for weakly nonlinear charging dynamics with u = 1, e = 0.05, 
5 = 0.1. The full numerical solution (solid) is compared with 
the approximate first-order expansion at the diffusion time 
scale, given by Eqs. 110811 . 110411 . and 11071 . Also, shown is 
the latter plus the zeroth-order inner approximation, Eq. I69I I. 
for the diffuse layers (dashed). 



of strength, —wq(pq)^ equal to the leading-order total 
salt adsorption. According to Eqs.lO and with 
jo(oo) = 0, this is given by 



wo(oo) = 4 sinh 



2 ,f-^(«) 



where 



/(C) = C + 2<5sinh(C/2) 



which reduces to 



Wq(co) — 4 sinh^ 



4' 



(109) 



(110) 



(111) 



in the absence of any compact layers (i5 = 0). 

This simple approximation describes two diffusion lay- 
ers created at the electrodes slowly invading the entire 
cell. At first, they have simple Gaussian profiles, 



5{x,t) - 1 



ewo(oo) 



,-(2;+l)V4t ^ g-(x-l)V4t 



(112) 



for f <C 1, which is qualitatively consistent with the nu- 
merical results in Fig. [7fc). To attempt a quantitative 
comparison, we also need t ^ 1 to use Eqs. H108|l and 
(|104|l . As shown in Fig. (Q, the approximation is rea- 
sonable for t = 3 with an error of roughly = 0.0025. 
The two diffusion layers eventually collide, and the con- 
centration slowly approaches a (reduced) constant value. 



c(a;, t) ~ 1 — ewoCoo) 



(113) 



for t ^ 1, as expected from the steady-state excess con- 
centration in the double layers. (This result may be 
checked by replacing the sum in Eq. I|104|) with an in- 
tegral in the limit t —* oo.) 



E. 



Bulk Concentration Polarization 



As mentioned above, the bulk charge density remains 
very small, p — 0{e^), even during double-layer relax- 
ation, but changes in neutral bulk concentration affect 
the potential at first order. Substituting the outer ex- 
pansions into Eq. H38() and collecting terms at 0(e), we 
have 



= 



d_ 

dx 



CO' 



dx 



Ci 



dx 



(114) 



This is easily integrated using co = 1 to obtain the first- 
order contribution to the bulk electric field, 



dx 



ii{t) - io{t) ci{x,t), 



(115) 



where the second term describes concentration polariza- 
tion, i.e. the departure from a harmonic potential, which 
would be predicted by Ohm's law. The first term is a 
uniform bulk field (or current) determined by first-order 
perturbation in double-layer charge. This follows from 
the matching condition, Eq. O, at first order: 



^91 -■ 



(116) 



where qi{t) is obtained by solving the inner problem at 
first order. 



F. Perturbations in Double-layer Structure 

Unfortunately, the first-order inner problem is difficult 
to solve analytically because the perturbed concentration 
profiles are no longer in thermal equilibrium during the 
initial charging phase. To see this, note that the time 
derivatives in Eqs. l(S7|) and contribute nonzero (but 
known) terms at first order. 



dcQ 


d 


dt 


dy 


dpo 


d 


dt 


dy 


d^i 




dy"^ 


= Pi 



dy 



dy 



dy 



dy 



(119) 



although one still solves a system of linear ordinary dif- 
ferential equations in y at each t, since CQ{y,t), p^{y,t), 
and (pQ^y^t) are known. 
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The general problem seems daunting, but some 
progress can be made at the scale of bulk diffusion, 
i = 0(1) or t = 0{e^^), where the leading-order concen- 
tration profiles remain in thermal equilibrium, without 
any explicit time dependence. This will give us some in- 
sight into secondary charge relaxation at the time scale 
of bulk diffusion. In this limit, the Equations (|117|l and 
()118|l can be integrated to obtain: 

~~d^ ^ Po(y,oo)— -f pi — (y,oo) (120) 

dpi . , d4>i d4)o, . , . 
-— ^ co{y,oo)—+ci — {y,oo) (121) 

after applying the usual van Dyke matching conditions. 
Substituting from Poisson's equation, p„ = — 3^0„/9y^, 
at orders n = 0, 1 into Eq. H120|) . integrating, and apply- 
ing matching again, we obtain: 

ci (y, t) ^ ^{y, oo) ^ (y, i) + ci (-1, t) (122) 

for t > and t — i/e 3> 1. From the previous section, we 
also have the leading-order inner concentration, 

co(2/,oo) = co(-l,<) + i-Bo(y,oo), (123) 

where co(— 1, f) = 1 is the leading-order outer concentra- 
tion, and 

Eo{y,oo) = — — y,oo) = 2smh , 124) 

dy 2 

is the leading-order inner electric field in steady-state. 
Finally, we substitute these expressions into Eq. H121|) 
and use Eq. (|118|l to obtain a master equation for the 

first-order inner electric field, Ei{y,t) = —^^{y,t), at 
the bulk-diffusion time scale: 



El + ciEq 



(125) 



This linear equation with a non-constant coefficient must 
be solved subject to the boundary conditions, Ei{oo,t) — 
and Ei{0,t) = —qi{t). The perturbation of the total 
charge, qi(t) is obtained by another integration of the 
field to get the first-order inner potential, while applying 
the Stern boundary condition. 

For our purposes here, it suffices to point out that 
the spatial profile of the first-order inner electric field in 
Eq. H125|) varies with the outer concentration, ci(— l,f), 
at the slow time scale of bulk diffusion. Notably, this 
can lead to a secondary relaxation of the total diffuse 
charge, in response to the evolution of the diffusion lay- 
ers. We observe this slow relaxation phase in our nu- 
merical solutions of the full equations, especially at large 
voltages. In particular, it is presumably associated with 
the non-monotonic charging profile for v = 4: shown in 
Fig. EIc). A detailed analysis of this interesting effect 
from the setup above would require solving the first-order 
inner problem numerically, so we leave it for future work. 



VIII. STRONGLY NONLINEAR DYNAMICS 
A. Steady State and the Dukhin Number 

We stress again that the asymptotic expansions de- 
rived above are valid in the limit of thin double layers, 
e — > 0, with the other two dimensionless parameters, v 
(applied voltage) and 5 (relative compact-layer capaci- 
tance), held constant. For any fixed e > 0, there is no 
guarantee that the approximation remains accurate as 
the other parameters are varied. Having just calculated 
the bulk concentration to first order in the regular ex- 
pansion, Eq. H53I) . wc can now check a posteriori under 
what conditions it remains a good approximation. 

A simple check involves the constant bulk concentra- 
tion, Eq. (|113|) . after the charging process is completed. 
The assumption that the first correction is much smaller 
than the leading term requires, as = ew(oo) <ti 1. Lin- 
earizing Eq. (|110(l for S 1, we can write this condition 
in a closed form: 



4 e sinh^ 



Ml + 6) 
Putting the units back, we have 



< 1 



MCo) = ^sinh^(^l«l 



(zeCo\ 



(126) 



(127) 



where (q w V/{1 + 6) is the steady-state zeta potential, 
long after the DC voltage is applied. 

The condition, asiCo) <C 1, for the validity of the 
steady-state asymptotic expansion is identical to that of 
small Dukhin number, Du(Co) ^ 1, from Eq. ^ in the 
limit of no electro-osmosis (m = 0), which may seem sur- 
prising since there is no surface conduction in our one- 
dimensional model problem. (Hence, we use the symbol 
as rather than Du.) The reason is that in both cases 
— Dukhin's problem of electrophoresis of highly charged 
particles in weak applied fields and ours of electrode 
screening in strong applied fields — the double layer ab- 
sorbs a significant amount of neutral salt from the bulk 
(reverse Donnan effect). 

Net charge adsorption relative to the point of zero 
charge is measured by the total zeta potential. 



Ctot Cgq ^" Cind 



(128) 



where Ceq is the uniform equilibrium zeta potential (re- 
ffecting the initial surface charge) and Qnd is the non- 
uniform induced zeta potential (resulting from diffuse- 
charge dynamics) . In Dukhin's problem, the former may 
be large, Du{Qq) > 1, but the latter is always small, 
Cind kT / ze, so that the charging dynamics is linearized 
(or ignored). In our model problem, the situation is re- 
versed: We assume Qq — (for simplicity), but we allow 
for a large applied voltage, Qnd ~ v/{l + 5) > kT/ze, in 
which case the dynamics is nonlinear. In both cases, the 
steady state is well described by weakly nonlinear asymp- 
totics as long as as(Ctot) = E)u(Ctot) ^ 1- When this 
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condition is violated, double-layer charging and surface 
conduction may cause significant changes in the steady- 
state bulk concentration. 



B. Breakdown of Weakly Nonlinear Asymptotics 

In general, weakly nonlinear dynamics breaks down 
at somewhat smaller voltages, where Qot > kT/e but 
Q^s(Ctot) = Du(Ctot) ^ 1, because neutral-salt adsorption 
causes a temporary, local depletion of bulk concentration 
exceeding that of the steady state, after diffusional re- 
laxation. In our model problem, the maximum change in 
bulk concentration occurs just outside the diffuse layers 
at a; = ±1, just after the initial charging process finishes 
at time scale, t = 1 or i = e. From Eq. (|112|l . we have 
the first two terms of the weakly nonlinear asymptotic 
expansion there: 



c(±l,e)^l 



wo(oo). 



(129) 



At that time, the newly created diffusion layers have 
spread to 0(-\/e) width, so the concentration is depleted 
locally by 0{e/y/e) = 0{y/e), which is much more than 
the uniform 0(e) depletion remaining after bulk diffu- 
sion. 

Therefore, in order for the time-dependent correction 
term to be uniformly smaller than the leading term, we 
need 



ad = J- t«o(oo) = 



< 1. 



The relevant dimensionlcss parameter. 




(130) 



(131) 



is larger than as{Ctot) (and the Dukhin number) by a 
factor of ^L/ttXd in the limit of thin double layers. For 
weakly nonlinear dynamics to hold, the applied voltage 
cannot greatly exceed the thermal voltage. 



V kT ^ L 

(tot ~ < log — , 

1 + ze Ad 



(132) 



even for very thin double layers. Ad <C L, due to the loga- 
rithm. In comparison, the applied voltage can be twice as 
large before the steady-state bulk concentration is signif- 
icantly affected (and surface conduction becomes impor- 
tant in higher dimensions). This could have interesting 
cons eque nces for induced-charge electrokinetic phenom- 
ena |12C| at moderate applied voltages where ad > I but 
a, = Du < 1. 



double-layer charging is coupled to bulk diffusion. As 
long as ad < 1, however, the bulk remains quasi- neutral 
at all times. This regime of strongly nonlinear dynamics 
is demonstrated by the numerical solution in Fig. 1101 for 
V — 4, e — 0.05, and 6 — 0.1, in which case ad — 0.545. In 
spite of the substantial 0(1) amount of charge transfered 
from one diffuse layer to the other, each retains almost 
exactly the same 0(e) width as at lower voltages, and 
bulk electro-neutrality remains an excellent approxima- 
tion for all times. The initial charging process up to t « 1 
creates a diffusion layer of neutral salt which relaxes into 
the bulk at the scale i « 1 (or t — i/e^ 20). 

In the strongly nonlinear regime, if ag is not too small, 
double-layer charging is slowed down so much by nonlin- 
earity that it continues to occur as the bulk diffusion 
layers evolve. One way to see this is that the effective 
RC time for the late stages of charging in Ea. (|85|l is 

t^{v) = Ci{v) « cosh I w 2sinh2 \~Yt ^^^^^ 

where we use the leading-order approximation of the dif- 
ferential capacitance, Eq. H88|l . for (5 ^ 1. In units of 
the bulk diffusion time, the nonlinear relaxation time is 
tc = etc — tts/2. 

To make analytical progress, one would consider the 
joint limits 

e ^ and u ^ oo with ad{v) > fixed (134) 

and expect the approximations to remain acceptable at 
somewhat larger voltages, as long as{v) < 1. Such anal- 
ysis is beyond the scope of this article, but at least we 
indicate how the leading order approximation would be 
calculated. (Going beyond leading order seems highly 
nontrivial.) 

At leading order in the bulk, we have the usual equa- 
tions for a neutral binary electrolyte (with equal ionic 
diffusivites) , 



dcp 

dt 



5% 



and 



d_ 

dx 



Co 



dx 



= 



(135) 



with p = O(e^). Integrating the second equation, we 
obtain a constant, uniform current density, joit), as be- 
fore, but the electric field is modified by concentration 
polarization. 



9(1)0 joit) 



dx co{x,t) 



(136) 



The effective boundary conditions come from asymptotic 
matching with the diffuse layers as before. 



C. Strongly Nonlinear Asymptotics 

When condition (|132|l is violated, electrochemical 
relaxation becomes much more complicated because 



dqo -. f.. , dwo dco 

only now the diffusive flux entering the diffuse layers 
(second equation) appears at leading order. The ionic 
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FIG. 10: Strongly nonlinear charging dynamics for v = 4 with e — 0.05 and S = 0.1. The potential (a), charge density (b), 
and concentration are shown in the half cell (top) and in the diffuse layer (bottom) for t = 0.5 (solid), 1 (dot), 2 (dash), 4 
(dot-dash), 8 (dot-dot-dot-dash), and 20 (long dash). 



concentrations retain Gouy-Chapman equilibrium pro- 
files modified quasi-statically by the evolving nearby bulk 
concentration: 

qo{t) = -2Vco(-l,i)sinh (^^^ (138) 

wait) - 4Vco(-l,t)sinh^ {^^^ ^^^^^ 

where 

Co{t) - qoit)S = ^o{t) = -V - 4(-l,t). (140) 

It seems exact solutions are not possible in terms of ele- 
mentary functions. The equations are "stiff", since they 
involve a short time scale, t — e, for the initial phases of 
charging, but at least the spatial boundary layers have 
been "integrated out" , which is convenient for numerical 
solutions. 



D. Space Charge at Very Large Voltages 

We close this section by noting some intriguing, new 
possibilities, further into the strongly nonlinear regime. 



At large voltages, such that 0;^ > 1, it seems a transient 
space charge layer should form since the bulk concentra- 
tion would be depleted almost completely near the dif- 
fuse layers by the initial charging process. In steady-state 
problems of Faradaic conduction, it is well known that 
double-layer structure is altered from its Gouy-Chapman 
equilibrium profile at a limiting current .132i | and may 
turn into an exten ded space charge layer above a lim- 
iting current |l33j . but here we see that similar effects 
may also occur temporarily with large time-dependent 
voltages, in the absence of any Faradaic processes (at 
blocking electrodes). At still larger voltages, such that 
Q!s > 1, double layer charging consumes most of the bulk 
concentration, presumably leaving the entire bulk region 
in a state of "space charge" . 

Such situations may seem quite exotic in macroscopic 
systems, where e = Xd/L is extremely small, but in mi- 
crosystems perhaps they could occur. The mathematical 
model neglects bulk reactions (e.g. leading to hydro- 
gen bubble formation), nonlinear dielectric properties, 
electro-convection, or other effects which may hinder the 
formation of space charge in real systems. Nevertheless, 
the rich nonlinear behavior of the model merits further 
mathematical study, as a challenging problem in time- 
dependent boundary-layer theory. 
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IX. BEYOND THE MODEL PROBLEM 

We conclude by discussing more general situations, 
which contain some new physics, absent in our simple 
model problem. For thin double layers, the same meth- 
ods of asymptotic analysis could be applied to derive ef- 
fective equations in which the double layers are incorpo- 
rated into boundary conditions, better suited for analyti- 
cal or numerical work. Here, we simply sketch the results 
and suggest some other model problems for further study. 

A. Two or More Dimensions 

In the weakly nonlinear regime, where < 1 for all 
times over all double layers, our analysis extends trivially 
to higher dimensions, as long as the surface curvature 
does not introduce another length scale much smaller 
than L. In that case, the double layers are locally "flat" , 
and the boundary-layer calculations remain unchanged. 
Following the same procedure, we find that the bulk con- 
centration is uniform at leading order, cq = 1, and the 
bulk potential, (j)(r,t), is a harmonic function, 

V^^o = 0, (141) 

subject to a (dimensionless) RC boundary condition at 
each electrode surface, 

— = 6 (00 - 0e) g-^ n • V<pQ, (142) 

where n is the unit normal pointing into the electrolyte 
and (l>e{^, t) is the local electrode potential relative to the 
solution. The latter is equal to the local applied voltage 
plus the equilibrium zeta potential: 

Mr,t) = Vir,t) + Ceg{r) (143) 

which accounts for any pre-existing double-layer charge 
(neglected in our calculations above). A Neumann 
boundary condition, n • V0Oj is imposed at any inert, 
non-polarizable surface, such as a channel side wall. 

Another complication in two or more dimensions is the 
possibility of electro-osmotic flow. The fluid velocity in 
the bulk usually satisfies the Stokes equations, which may 
be unsteady for high-frequency forcing. In the weakly 
nonlinear regime, the classical Helmholtz-Smoluchowski 
formula gives the fluid slip in terms of the local zeta po- 
tential and tangential bulk electric field |23, 12M l29| . 

Equations H141|) and H142|) model the electrolyte as 
a bulk Ohmic resistor with a capacitor skin at elec- 
trode interfaces. The linearized version of these equations 
(with C — constant) has been studied extensively, e.g. in 
the context of metallic colloids |lO0l Il04| . AC electro- 
osmosis AC pumping Wl^d^er phenomena 
of induced-charge electro-osmosis |l20l Il2lj . The nonlin- 
ear version, however, has apparently not been analyzed, 
even though it may have relevance for experiments, in 



which the condition, v <^ 1 {V kT/ze), is routinely 
violated. 

More significant modifications arise at leading order in 
the strongly nonlinear regime (or at higher order in the 
weakly nonlinear regime). Ohm's law breaks down due 
to concentration gradients, as the double layers absorb 
a significant amount of neutral salt from the bulk. In 
two or more dimensions, the dimensionless leading-order 
equations for co(r, t) and ^o(r, t) in section rVlIII take the 
form, 

^^V^co, V.(coV0) = O, (144) 
ot 

where we scale time to the bulk diffusion time. This 
assumes ad < 1 so that no transient space charge layers 
form. 

The effective boundary conditions still involve the 
small parameter, e, as in one dimension, since the natural 
scale is the RC charging time, but there are some new 
terms in higher dimensions: 

= n.(coV,^o)-DuV,- J, (145) 
= n-Vco^T)uVs-{DsVsWo) (146) 

ot 

The last term in Eq. (|145|l is the surface divergence of 
the (leading-order) dimensionless tangential current, J^, 
in the diffuse layer; the size of this term compared to the 
normal current is governed by a Dukhin number, based 
on the largest expected total zeta potential. Similarly, 
the last term in Eq. H146|) is the surface divergence of the 
tangential diffusive fiux in the diffuse layer, where Dg is 
a dimensionless surface diffusivity; again, this term is of 
order Du smaller than the normal diffusive fiux. 

Formulae for Jg and Dg can be derived systematically 
using the matched asymptotic expansions, which is be- 
yond the scope of this paper. The classical r esult s of 
Bikerman [lid Ill7| and Deryagin and Dukhin |l37j are 
available for the case of weak applied voltages {Qnd ^ 
kT/e) and large equilibrium surface charges {Qq > kT/e, 
Du(Ceg) ~ 1), and many Russian authors have stud- 
ied electrokinetic phenomena in this regime [H, l33| . 
The case of strongly nonlinear dynamics {(ind > kT/e, 
«d(Cmd) ~ 1), however, should be revisited in more de- 
tail to see if any changes arise for strong, time-dependent 
applied voltages. We suggest as a basic open question an- 
alyzing the electrochemical response of a metal cylinder 
or sphere in a strong, suddenly applied, uniform back- 
ground DC field. 

Another interesting issue is the stability of our one- 
dimensional solution. One should consider small space- 
dependent perturbations of the solution at various large 
voltages, in both the weakly nonlinear and strongly non- 
linear regimes. The general transient analysis in two or 
more dimensions with the same equations and boundary 
conditions presents an interesting challenge. 
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B. General Electrolytes and Faradaic Reactions 

Even in one dimension, it would be interesting to ex- 
tend our analysis to more general situations involving 
asymmetric or multicomponent electrolytes, which un- 
dergo Faradaic processes at electrode surfaces. Restor- 
ing dimensions, the bulk electrolyte is described by the N 
ionic concentrations, Ci, i = 1, 2, . . . , A'^, satisfying mass 
conservation, 



where is the flux density due to diffusion and electro- 
migration, 

= -AVC, - mz,eaV^ (148) 

as in Eq. For thin double layers, at leading order 

the bulk remains neutral (as long as < 1 to avoid 
space charge formation), so the potential is determined 
implicitly by the condition of clectroneutrality, 

N 
i=l 

These are the standard equations of bulk electrochem- 
istry [23, but interesting physical effects are contained 
in the effective boundary conditions. 

Generalizing the total surface charge density q and ex- 
cess surface concentration w, we define to be the sur- 
face concentration of species i absorbed in the diffuse 
layer. To be precise, it is the integral of the leading- 
order excess concentration relative to the bulk over the 
inner coordinate, as in Eqs. H71|l and 1)93(1 . For example, 
q = ze(F-|_ — r_)/2 and w = (F-|_-|-F_)/2 for a symmetric 
binary electrolyte. 

Following the procedure above, the boundary condi- 
tions on the leading-order bulk approximation are of the 



form: 




= n • F, + V, • F,, + R, (150) 



where Fsi{Ci, $) is th e surface flux density of species i in 
the double layer |l37j and Ri{{Ci}, $) is the reaction-rate 
density for any Faradaic processes consuming (or produc- 
ing) species i at the surface. The usual assumption for 
Ri involves Arrhenius kinetics, as in the Butler- Volmer 
equation, but the Frumkin correction for concentration 
variations across the diffuse layer must be taken into ac- 
count mEi. 

The general system of nonlinear equations is challeng- 
ing to solve, even numerically, due to multiple length and 
time scales. Boundary-layer theory provides only a par- 
tial simplification by integrating out the smallest length 
scale. As described in section |nl various special cases 
of the effective equations have been considered in the 
literature, but much remains to be done, especially for 
strongly nonlinear dynamics in large applied voltages. In 
micro-electrochemical or biological systems, this regime 
is easily reached, so it merits additional mathematical 
study and comparison with experimental data, in part 
to test the applicability of the Nernst-Planck equations 
in micro-systems. Another interesting aspect is the cou- 
pling of electrochemical dynamics to fluid flow, which is 
finding new applications in microfluidic devices. 



Acknowledgments 

This work was supported in part by the MRSEC Pro- 
gram of the National Science Foundation under award 
number DMR 02-13282 (MZB), by ESPCI through the 
Paris Sciences Chair (MZB), and by an ACI from the 
Ministere de la Recherche (AA). We are grateful to K. 
Chu and the MIT librarians for help in obtaining articles 
from the "ancient" and Russian literatures. 



[1] A. Ramos, H. Morgan, N. G. Green, and A. Castellanos, 

J. Phys. D: Appl. Phys. 31, 2338(1998). 
[2] A. Ramos, H. Morgan, N. G. Green, A, Castellanos, J. 

Colloid Interface Sci. 217, 420(1999). 
[3] N. G. Green, A. Ramos, A. Gonzalez, H. Morgan, and 

A. Castellanos, Phys. Rev. E 61, 4011 (2000). 
[4] A. Gonzalez. A. Ramos, N. G. Green, A. Castellanos, 

and H. Morgan, Phys. Rev. E 61, 4019(2000). 
[5] N. G. Green, A. Ramos, A. Gonzalez, H. Morgan, and 

A. Castellanos, Phys. Rev. E 66, 026305 (2002). 
[6] A. Ajdari, Phys. Rev. E 61, R45(2000). 
[7] A.B.D. Brown, C.G Smith, A.R. Rennie, Phys.Rev.E, 

63, 016305 (2001). 
[8] V. Studer, A. Pepin, Y. Chen, and A. Ajdari, Microelec- 
tronic Eng. 61, 915 (2002). 
[9] M. Mpholo, C.G. Smith, and A.B.D. Brown, Sensors and 

Actuators B 92, 262 (2003). 



[10] A. Ramos, A. Gonzalez, A. Castellanos, N. G. Green, 
and H. Morgan, Phys. Rev. E 67, 056302 (2003). 

[11] F. Nadal, F. Argoul, P. Kestener, B. Pouligny, C. Ybert, 
and A. Ajdari, Eur. Phys. J. E 9, 387 (2002). 

[12] S. Yeh, M. Seul, and B. Shraiman, Nature 386, 57 
(1997). 

[13] M. Trau, D. A. Saville, and I. A. Aksay, Langmuir 13, 
6375 (1997). 

[14] C. Faure, N. Decoster, and F. Argoul, Eur. Phys. J. B 5, 
87 (1998). 

[15] N. G. Green, A. Ramos, and H. Morgan, J. Phys. D.: 

Appl. Phys. 33, 632 (2000). 
[16] F. Nadal, F. Argoul, P. Hanusse, B. Pouligny, and A. 

Ajdari, Phys. Rev. E 65, 061409 (2002). 
[17] C. Marquet, A. Buguin, L. Talini, and P. Silberzan, Phys. 

Rev. Lett. 88, 168301 (2002). 
[18] W. D. Ristenpart, I. A. Aksay, and D. A. Saville, Phys. 



25 



Rev. Lett. 90, 12803 (2003). 

W. Helfrich, Z. Naturforsch C 29, 182 (1974). 

M. D. Mitov, P. Meleard, M. Winterhalter, M. I. An- 

gelova, and P. Bothorel, Phys. Rev. E 48, 628 (1993). 

R. Pethig, Critical Rev. Biotechnol. 16, 331(1996). 

L. Journiaux and J. P. Pozzi, J. Geophys. Res. 102, 

15335 (1997). 

P. Reppert, F. D. Morgan, D. P. Lesmes, and L. Jouni- 
aux, J. Colloid Interface Sci. 234, 194 (2001). 
P. M. Reppert and F. D. Morgan, J. Colloid Interface 
Sci. 254, 373 (2002). 

A. J. Bard and L. R. Faulkner, Electrochemical Methods: 
Fundamentals and Applications (Wiley, New York, 1980). 
J. Newman, Electrochemical Systems (Prentice-Hall, En- 
glewood Cliffs, NJ, 1991). 

R. J. Hunter, Foundations of Colloidal Science (Oxford, 
second edition, 2000). 

W. B. Russel, D. Saville, and W. R. Schowalter, Colloidal 
dispersions (Cambridge University Press, 1989). 
J. Lyklema, Fundamentals of Interface and Colloid Sci- 
ence, Vol. 2 (Academic Press, New York, 1995). 
J. D. Ferry, J. Chem. Phys. 16, 737 (1948). 
J. R. Macdonald, Trans. Faraday Soc. 67, 943 (1970). 
A. A. Kornyshev, and M. A. Vorotyntsev, Electrochimica 
Acta 26, 303 (1981). 

S. S. Dukhin and V. N. Shilov, Adv. Colloid Inteface Sci. 
13, 153 (1980). 

S. S. Dukhin, Adv. Colloid Interface Sci. 44, 1 (1993). 
M. Scott, K. V. I. S. Kaler, and R. Paul, J. Colloid In- 
terface Sci. 238, 449 (2001). 

A. Ramos, A. Gonzalez, N. G. Green, H. Morgan, and 
A.Castellanos, J. Colloid Interface Sci. 243, 265 (2001). 
J. R. Macdonald, Electrochim. Acta 35, 1483 (1990). 
L. A. Geddes, Ann. Biomedical Eng. 25, 1 (1997). 

C. M. Bender and S. A. Orszag, Advanced Mathemat- 
ical Methods for Scientists and Engineers: Asymptotic 
Methods and Perturbation Theory (Springer, 1999). 

E. J. Hinch, Perturbation Methods (Cambridge Univer- 
sity Press, 1991). 

J. Kevorkian and J. D. Cole, Multiple Scale and Singular 
Perturbation Methods (Springer, New York, 1996). 

F. Kohlrausch, Pogg. Ann. 148, 143 (1873). 

H. L. F. von Helmholtz, Ann. Physik 89, 211 (1853). 
H. L. F. von Helmhohz, Ann. Phys. Chem. 7, 337 (1879). 

G. Gouy, Ann. Chim. Phys. 29, 145 (1903). 

E. Warburg, Ann. Phys. Chim. 67, 493 (1899). 

E. Warburg, Ann. Phys. 6, 125 (1901). 

F. Kriiger, Z. Phys. Chem. 45, 1 (1903). 

J. E. B. Randies, Disc. Faraday Soc. 1, 11 (1947). 

G. Gouy, Compt. Rend. 149, 654 (1909). 
G. Gouy, J. Phys. (4) 9, 45 (1910). 

A. Einstein, Ann. Physik 19, 371 (1905). 

P. Debye and E. Hiickel, Physik. Z. 24, 185 (1923). 

P. Debye and E. Hiickel, Physik. Z. 24, 305 (1923). 

D. L. Chapman, Phil. Mag. 25, 475 (1913). 

P. Delahay, Double Layer and Electrode Kinetics (John 
Wiley, New York, 1965). 

J. O'M. Bockris and A. K. N. Reddy, Modern Electro- 
chemistry (Plenum, New York, 1970). 
M. Sluyters-Rehbach and J. H. Sluyters, in Electroana- 
lytical Chemistry 4, 1-128, ed. by A. J. Bard (Marcel 
Dekker, New York, 1970). 
R. Parsons, Chem. Rev. 90, 813 (1990). 
O. Stern, Z. Elecktrochem. 30, 508 (1924). 



[61] M. A. Vorsina and A. N. Frumkin, Compt. Rend. Acad. 

Sci. U. R. S. S. 24, 918 (1939). 
[62] D. C. Grahame, Chem. Rev. 41, 441 (1947). 
[63] J. R. Macdonald, J. Chem. Phys. 22, 1857 (1954). 
[64] D. C. Grahame, J. Chem. Phys. 18, 903 (1950). 
[65] A. D. HoUingsworth and D. A. Saville, J. Colloid and 

Interface Sci. 257, 65 (2003). 
[66] I. L. Cooper and J. A. Harrison, J. Electroanal. Chem. 

86, 425 (1978). 
[67] W. Nernst, Z. Physik. Chem. 2, 613 (1888). 
[68] W. Nernst, Z. Physik. Chem. 4, 129 (1889). 
[69] M. Planck, Ann. Phys. Chem. 39, 161 (1890). 
[70] G. Jaffe and C. Z. LeMay, J. Chem. Phys. 21, 920 (1953). 
[71] A. A. Kornyshev, and M. A. Vorotyntsev, Phys. Stat. 

Sol. (a) 39, 573 (1977). 
[72] A. A. Kornyshev, and M. A. Vorotyntsev, Electrochimica 

Acta 23, 267 (1978). 
[73] R. P. Buck, J. Electoanal. Chem. 23, 219 (1969). 
[74] P. Debye and H. Falkenhagen, Physik. Z. 29, 121 (1928). 
[75] H. Falkenhagen, Electrolytes (Oxford University Press, 

1934). 

[76] G. Jaffe, Ann. Physik 16, 249 (1933). 
[77] G. Jaffe, Phys. Rev. 85, 354 (1952). 
[78] H.-C. Chang and G. Jaffe, J. Chem. Phys. 20, 1071 
(1952). 

[79] J. R. Macdonald, Phys. Rev. 92, 4 (1953). 

[80] R. Friauf, J. Chem. Phys. 22, 1329 (1954). 

[81] C. G. J. Baker and E. R. Buckle, Trans. Faraday Soc. 
64, 469 (1968). 

[82] J. R. Macdonald, J. Chem. Phys. 61, 3977 (1974). 

[83] J. R. Macdonald, J. Electroanal. Chem. 53, 1 (1974). 

[84] J. H. Beaumont and P. W. M. Jacobs, J. Phys. Chem. 
solids 28, 657 (1967). 

[85] E. M. Itskovich, A. A. Kornyshev, and M. A. Vorotynt- 
sev, Phys. Stat. Sol. (a) 39, 229 (1977). 

[86] J. Gunning, D. Y. C. Chan, and L. R. White, J. Colloid 
interface Sci. 170, 522 (1995). 

[87] M. Scott, R. Paul, and K. V. I. S. Kaler, J. Colloid In- 
terface Sci. 320, 377 (2000). 

[88] M. Scott, R. Paul, and K. V. I. S. Kaler, J. Colloid In- 
terface Sci. 320, 388 (2000). 

[89] E. H. B. DeLacey and L. R. White, J. Chem. Soc. Fara- 
day Trans. 2 78, 457 (1982). 

[90] V. G. Levich, Doklady Akad. Nauk SSSR 67, 309 (1949). 

[91] A. Frumkin, Z. Electrochem. 59, 807 (1955). 

[92] J. R. Macdonald, J. Electroanal. Chem. 70, 17 (1976). 

[93] J. R. Macdonald and P. W. M. Jacobs, J. Phys. Chem. 
Solids 37, 1117 (1976). 

[94] D. R. Franceschetti and J. R. Macdonald, J. Electroanal. 
Chem. 82, 271 (1977). 

[95] A. Bonnefont, F. Argoul, and M. Z. Bazant, J. Elec- 
troanal. Chem. 500, 52 (2001). For a preprint with more 
results from Ref. [gg, see cond-mat/0006104. 

[96] A. Bonnefont, Local methods of characterizing the dy- 
namics of transport and reactivity of an interface, Ph.D. 
thesis, Universite de Bordeaux (2000). 

[97] T. Brumleve and R. Buck, J. Electroanal. Chem 90, 1 
(1978). 

[98] W. D. Murphy, J. A. Manzanares, S. Male, and H. Reiss, 

J. Phys. Chem. 96, 9983 (1992). 
[99] H. Cohen and J. W. Cooley, Biophys. J. 5, 145 (1965). 
[100] V. A. Murtsovkin, KoUoidn. Zh. 58, 358 (1996). 
[101] S. Dukhin, in Encyclopedia of Surface and Colloid 

Sience, ed. by M. Dekker, 3677 (2002). 



26 



[102; 
[103; 

[104 

[105; 
[106; 
[107; 
[los; 

[109 
[110' 

[111 

[112' 
[113^ 

[114; 

[115" 
[116 
[117^ 
[118 
[119^ 
[120^ 

[121 

[122; 

[123' 

[124 
[125" 
[126^ 

[127 
[128 
[129 

[130 

[131 
[132; 
[133; 
[134; 
[135; 

[136; 

[137 
[138 



V. G. Lcvich, Physicochemical Hydrodynamics 
(Prentice-Hall, Englewood Cliffs, NJ, 1962). 

I. N. Simonov and V. N. Shilov, Kolloidn. Zh. 35, 381 
(1973). 

I. N. Simonov and V. N. Shilov, Kolloidn. Zh. 39, 878 
(1977). 

N. I. Gamayunov, V. A. Murtsovkin, and A. S. Dukhin, 
Kolloidn. Zh. 48, 233 (1986). 

I. N. Simonov and S. S. Dukhin, Kolloidn. Zh. 35, 191 
(1973). 

T. S. Simonova and S. S. Dukhin, Kolloidn. Zh. 38, 86 
(1976). 

V. N. Shilov and T. S. Simonova, Kolloidn. Zh. 43, 114 

(1981). 

T. S. Simonova, V. N. Shilov, and O. A. Shramko, Kol- 
loidn. Zh. 63, 114 (2001). 

B. V. Deryagin, S. S. Dukhin, and V. N. Shilov, Adv. 
Colloid. Interface Sci. 13, 141 (1980). 

S. S. Dukhin and B. V. Derjaguin, in Surface and Col- 
loidal Science, od. by E. Matijevic and F. R. Eirich, Vol. 

7 (Wiley, New York, 1974). 

R. W. O'Brien, J. Colloid Interface Sci. 92, 204 (1983). 

E. J. Hinch, J. D. Sherwood, W. C. Chew, and P. N. 
Sen, J. Chem. Soc, Faraday Trans 2 80, 535 (1984). 

D. C. Prieve, J. L. Anderson, J. P. Ebel, and M. E. 
Lowell, J. Fluid Mech. 148, 247 (1984). 

J. L. Anderson, Ann. Rev. Fluid Mech. (1989). 

J. J. Bikcrman, Z. Physik. Chcm. A163, 378 (1933). 

J. J. Bikcrman, KoUoid Zcitschrift 72, 100 (1935). 

J. J. Bikerman, Trans. Faraday Soc. 36, 154 (1940). 

J. T. G. Overbeek, KoUoidchem. Beih. 59, 287 (1943). 

M. Z. Bazant and T. M. Squires, preprint, 
http://arXiv.org/abs/physics/0306100. 

T. M. Squires and M. Z. Bazant, preprint, 
http://arXiv.org/abs/physics/0304090. 

B. M. Grafov and A. A. Chernenko, Dokl. Akad. Nauk. 
S.S.S.R. 146, 135 (1962). 

A. A. Chernenko, Dokl. Akad. Nauk SSSR, 153(5), 1129 
(1963). 

J. Newman, Trans. Faraday Soc. 61, 2229 (1966). 

A. D. Macgillivray, J. Chem. Phys. 48, 2903 (1968). 

I. Rubinstein, Electro- Diffusion of Ions (SIAM, 
Philadelphia, 1990). 

J. Henry and B. Louro, Nonlinear Anal. 13, 787 (1989). 

B. Louro, Portugaliae Mathcmatica 48, 179 (1991). 

J. Henry and B. Louro, Asymptotic Anal. 10, 297 
(1995). 

J.-H. Park and J. W. Jerome, SIAM J. Appl. MAth. 57, 
609 (1997). 

V. Barcilon, D.-P. Chen, R. S. Eisenbcrg, and J. W. 
Jerome, SIAM J. Appl. Math. 57, 631 (1997). 

W. Smyrl and J. Newman, Trans. Faraday Soc. 63, 207 
(1967). 

I. Rubinstein and L. Shtilman, J. Chem. Soc. Faraday 
Trans. II 75, 231 (1979). 

S. S. Dukhin and V. N. Shilov, Kolloidn. Zh. 31, 706 

(1969). 

S. S. Dukhin and N. M. Semenikhin, Kolloidn. Zh. 32, 

360 (1970). 

R. W. O'Brien and L. R. White, J. Chem. Soc. Faraday 
Trans. II 74, 1607 (1978). 

B. V. Deryagin and S. S. Dukhin, Kolloidn. Zh. 31, 350 

(1969). 

F. G. Donnan and A. B. Harris, J. Chem. Soc. 99, 1554 



(1911). 

[139] F. G. Donnan, Z. Eloktrochcm. 17, 572 (1911). 
[140] F. G. Donnan, Chcm. Rev. 1, 73 (1924). 
[141] J. Lyklema, Fundamentals of Interface and Colloid Sci- 
ence, Vol. 1 (Academic Press, New York, 1995). 



